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Abstract 

Game theory provides a general mathematical background to study the effect of pair interactions and 
evolutionary rules on the macroscopic behavior of multi-player games where players with a finite number 
of strategies may represent a wide scale of biological objects, human individuals, or even their associations. 
In these systems the interactions are characterized by matrices that can be decomposed into elementary 
matrices (games) and classified into four types. The concept of decomposition helps the identification of 
potential games and also the evaluation of the potential that plays a crucial role in the determination of 
the preferred Nash equilibrium, and defines the Boltzmann distribution towards which these systems evolve 
for suitable types of dynamical rules. This survey draws parallel between the potential games and the 
kinetic Ising type models which are investigated for a wide scale of connectivity structures. We discuss 
briefly the applicability of the tools and concepts of statistical physics and thermodynamics. Additionally 
the general features of ordering phenomena, phase transitions and slow relaxations are outlined and applied 
to evolutionary games. The discussion extends to games with three or more strategies. Finally we discuss 
what happens when the system is weakly driven out of the ’’equilibrium state” by adding non-potential 
components representing games of cyclic dominance. 
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1. INTRODUCTION 


Most games can be considered as simplified real-life situations in which the number of players as well as 
their options is limited. The simplest case is represented by two-player non-cooperative games where the 
intelligent and selfish players wish to maximize their own payoff that is quantified by a real number and 
depending on their simultaneous choices while their utilities are quantified by two payoff matrices [lj. 

The original concept and application of payoff matrices have been extended basically within the frame¬ 
work of evolutionary game theory U i i S i 0, i i 0- In biological systems the payoff matrices 
characterize the resultant fitness (more precisely, the capability to create offspring) for the interacting pair 
of species (representing the strategies for this terminology) and serve as a fundamental quantity in the 
population dynamics developed in the spirit of the Darwinian selection PEEMl. In human systems, 
however, the imitation of the more successful player can control the time-dependent frequency of players 
following a given strategy, in close analogy to biological systems eh, 3 , ii m. 

The first systematic investigations of evolutionary games were performed in well-mixed population of 
players for a wide scale of games and strategies. A progressively expanding research area was initiated by 
Nowak and May |20[ who introduced a model where the players are distributed on a square lattice and 
their incomes come from one-shot games with their neighbors. In this cellular automaton type model the 
players changed their strategy simultaneously (in discrete time steps) by imitating the strategy of the most 
successful neighbor. The numerical investigations of these deterministic models 21] have demonstrated the 
advantage of spatial structures for short-range interactions in the maintenance of cooperation among selfish 
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players for the prisoner’s dilemma. Since that time numerous followers have clarified relevant effects and 
phenomena supporting the maintenance of cooperative behaviors in real-life situations described by different 
versions of social dilemmas as surveyed by Nowak [gj, Sigmund Szabo and Fath [22], Allen and Nowak 


The systematic investigations involve systems where the authors studied the effects of different sets of 
strategies, the connectivity structures described by lattices or graphs, and a wide scale of possible evolu¬ 
tionary rules including co-evolutionary games where all the relevant ingredients of the model are allowed to 
evolve together with the strategy distribution (24j. Up to now most of the relevant two-person games are 
studied within the framework of evolutionary games. New features can be inv estig ated via versions of group 
interactions that cannot be decomposed into the sum of pair interactions 
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Within this wide scale of systems less efforts are focused on the potential games [27], which are related 
intimately to the statistical physics 128]. For a more detailed exposition of potential games see the relevant 
chapter in the book by Sandholm [101 ]. The decomposition of the matrix games into the linear combination 
of elementary games was discussed previously by Candogan et al. Q, Hwang and Rey-Bellet [3C$. The 
latter approach allows us to identify different types of pair interactions and also to determine the existence 
of potential as well as to evaluate it if it exists [3l|. Very recently Cheng 32] has suggested another method 
for the evaluation of the potential. The applications and perspectives of potential games within the economy 
are reviewed by Mallozzi 33]. 

In this review we survey briefly the fundamental concepts and background material necessary for under¬ 
standing the recent and future developments in the research of evolutionary potential games. Due to the 
wide scale of research fields involved and the limited length of this article we cannot give a comprehensive 
picture. Instead of mathematical rigorousness we use a concise style with an effort to give a wide picture of 
the set of phenomena and relationships that are already studied in the fields of statistical physics. 

The structure of this review is as follows. The next four Sections summarize the basic concepts and 
methods of traditional and evolutionary game theory including the general features of potential games and 
the way of the evaluation of potential. Subsequently we discuss the intimate relationship between multi¬ 
agent evolutionary potential games and Ising type models. Afterwards we briefly survey the phenomena (e.g., 
order-disorder phase transitions) explored by the application of the Ising type models. In Sec. [5]we describe 
briefly the ordering processes characterizing the ways how the systems tend towards the final stationary 
ordered structure at low noise levels. Section El is addressed to discuss what happens when the system 
behavior deviates from the thermodynamical equilibrium due to the presence of non-potential components 
representing cyclic dominance in the payoff matrices. Finally we outline some challenging questions for the 
continuation of research in the near future. 


2. BRIEF SURVEY OF GAMES 

First we give a concise and tutorial description of the basic definitions and concepts of normal games we 
use throughout this work. The reader can find further details of the traditional game theory in standard 
textbooks like Fudenberg and Tirole Q, Gibbons [35j, Hofbauer and Sigmund @, Weibull [3a], Samuelson 
f37| . Gintis Cressman Q, Sandholm [10], Sigmund [o] that provide a more general outline of a wider 
scale of games and illustrate their application in biology, economics, and behavioral sciences. 

2.1. Players, strategies, payoffs, and potential 

The normal game is an abstract formulation of a decision situation where each of N players must 
choose simultaneously (and independent of each other) one of their own possible options to receive payoffs 
dependent on the choices of all of them. Each intelligent and selfish (rational) player wishes to maximize 
her own payoff with the assumption that the co-players are also intelligent. In other words, the players 
are capable of deducing the best possible way of playing the game and they can handle the problem of the 
common knowledge implying the fact that all players know that the others wish to maximize their own 
utility, and all players know that all the others know that all of them wish to maximize their own payoff, 
etc. 
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In multi-player games the players will be denoted by x = 1,... ,7V. Each player x can choose a pure 
strategy s x from her set of strategies indexed by an integer, s x £ {1, 2 ,..., n x } where n x denotes the number 
of strategies of player x. The choices of all the players are defined by the strategy profile s = (si,..., Sn ) 
determining their payoff. In this notation each player x receives a payoff u x ( s) = u x (si ,..., sn) quantifying 
her utility. The above strategies are pure strategies. In many games (e.g., poker, matching pennies, or 
rock-scissors-paper game) the players can also play mixed strategies where they select one of their pure 
strategies with a given probability in each decision instance. Henceforth our discussion is limited to the 
cases where the players are constrained to use pure strategies. 

In the above games each of the TV players wishes to optimize her own payoff regardless of the others, 
which is complicated by the fact that, in general, none of the players can achieve her maximum payoff 
because of the counter-interest of the others. Despite it, the players can be satisfied if they choose a strategy 
profile s* = (s*,..., s* N ), called Nash equilibrium, that satisfies the following conditions: 


Vz,Vs x / 4 : u x (s* x ,s *_ x ) > u x {s x ,s*_ x ). 


(1) 


where S- x = (s i,..., s x _i, s s +i,..., sn) denotes the strategy profile of the co-players. If the inequality is 
strict then s* is called a strict Nash equilibrium. For the Nash equilibrium the players do not have unilateral 
incentive to choose another strategy. In other words, each player is satisfied with her own choice, because 
there is no way to receive a higher payoff for the given choices of the others. Nash’s theorem 38), |39[ states 
that in normal games there exists at least one Nash equilibrium, possibly involving mixed strategies. 

For most of the cases the game has more than one Nash equilibrium. In these situations we need 
additional criteria for the strategy selection. For example we can choose those Nash equilibria that provides 
the higher sum of the individual utilities. Another way of selection, based on the consideration of risk 
dominance, was introduced by Harsanyi and Selten ^0i|. The latter method suggests choosing the strategy 
that yields higher expected payoff against an opponent playing a mixed strategy when all feasible strategies 
are chosen with equal probability. These additional criteria may give contradictory advices for some games. 

Some of the games have attracted huge attention due to their curiosity and the high frequency we face 
them day by day. The best known example is the prisoner’s dilemma representing a real-life situation when 
the two players have only two options, called cooperation and defection, and the game has only one Nash 
equilibrium, the mutual defection, while the mutual cooperation would result in a higher payoff for both 
players (hence the dilemma). The analysis of similar situations attracted pro gres sive activity in different 
fields of science, including biology, economics, social sciences, and physics [8l. Ifli. 12211 - 

In normal games each player has her own utility function u x ( s) to consider. For the potential games we 
can introduce a single function F(s), called potential, that involves the strategic incentives of all players 
in the following way suggested by Monderer and Shapley [27j. Namely, the difference of the potential 
F(s) = V(si, S 2 , ■ ■ ■, sn), when player x modifies her strategy (from s x to s' x ), is equal to the utility 
function difference for the given player, that is, 


u x (s' x , s- x ) - u x (s x ; s - x ) = V{s' x , s - x ) - V(s x ;S- x ) 


( 2 ) 


for \/ x, Vs x , s' x , and V S- x . This means that the potential is constructed from the payoff variations of those 
players changing their strategies. Evidently, the existence of a potential for a normal game is strongly 
limited by the above conditions. In other words, Eq. m gives strict constraints for the possible values 
of payoffs. The difficulties can be well illustrated by introducing a dynamical graph 4l| where each node 
represents a strategy profile s and the edges denote transitions between two strategy profiles when only a 
single player modifies her strategy. The sum of the given individual payoff differences is zero along all loops 
of this dynamical graph for potential games. 

Notice that if only one player modifies her strategy along a loop, e.g., s x (l) —► s x (2 ) s x (k ) —► 

s^l)) (while S- x is quenched) then the sum of the payoff variation: 


k 

y, [u x (s x (i + 1 ); s_ 2 ,) - U x (s x (i); s_ x )] = 0 (3) 

2=1 
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because in this sum all payoff components appear twice with opposite signs as s x (k + 1) equals s x (l). 
Monderer and Shapley 27] have proved that a normal game admits a potential if and only if over any four- 
edge loops of the dynamical graph the sum of the changes in the deviator’s payoffs equals zero. Analogously 
to the potential energy in physical systems the potential is unique also here up to addition of a constant. 
In the literature of game theory one can find a wide range of examples satisfying the conditions of potential 
games from genetic competition 0,1 to congestion game on network 0 , FTaL 0 ] or oligopoly games 
471. 

The existence of the potential simplifies the analysis. For example, the strategy profile where V(s) 
achieves its maximum is a preferred (pure) Nash equilibrium that plays a distinguished role in evolutionary 
game theory like the ground state in physical systems. 

If the edges of the above mentioned dynamical graph are directed (called flow graph henceforth) by 
pointing towards the higher potential value, then the nodes (strategy profiles) with only incoming edges are 
Nash equilibria. Similar flow graphs can be derived for the so-called ordinal potential games where a weaker 
condition should be satisfied. Namely, if 


u x {s' x \s- x ) - u x (s x ; s- x ) > 0 


( 4 ) 


then 


V{s' x -,s- x ) -V{s x \s- x ) > 0. 


( 5 ) 


Evidently, the existence of a potential prohibits the presence of directed loops here. The advantage of the 
latter feature will be demonstrated later. In the literature of game theory one can find other classes of 
potential games, that we will not discuss here. For example, besides the exact and ordinal potential games 
mentioned above, Monderer and Shapley 27] have introduced the weighted and generalized ordinal potential 
games; Voorneveld 0 has studied the best-response potential games; and Morris and Ui 0 investigated 
the robust sets of equilibria for the generalized potentials. 

The most relevant application of potential games was discovered by Blume 00, who proved that for 
a certain set of evolutionary rules the system evolves into a Boltzmann-Gibbs ensemble and the stationary 
state can be well investigated by the tools of statistical physics. This feature justifies the importance of the 
preferred Nash equilibrium. 

Henceforth our analysis will be focused on those multi-player games where the players are allowed to use 
only pure strategies and the interaction is composed of two-player games surveyed in the following section. 


2.2. Two-player games 

For a two-player game with players x and y the payoffs are defined by payoff tables and we can apply the 
terminology of matrices. In many cases the strategy labels i and j (i = 1,2,..., n x and j = 1,2,..., n y ) are 
sufficient to identify the pure strategies selected by players x and y. For the expression of payoffs, however, 
it is convenient to denote the strategies of player x by (n^-dimensional) unit vectors, as 


®:r — — 


I 1 ) 
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, 83,2 = 

1 
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n 

0 

Vi / 


( 6 ) 


and similar expressions define the strategies of player y. For this notation the payoffs of player x and y can 
be expressed by the following matrix products as 


^(^S,, Sy ) - S x ' AS y , 

^y(^3,, ^y) — Sy • BS„ — Sa, • B S y , 


( 7 ) 


where the matrix elements Aij and Bij define the payoffs of players x and y (we used the relation Z?7 = Bji). 
Here we have to mention that the above expressions define the average payoffs for mixed strategies, e.g., 
s x = '}2 i piS X i, when player x uses her zth strategy with a probability /?* (JA pi = 1). The latter notation 
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is useful for games having mixed Nash equilibrium. If the players are different (like in male-female, young- 
old, buyer-seller, sender-receiver interactions) then the so-called asymmetric two-player normal games are 
specified by two payoff matrices (A and B) and it is customary to use its bi-matrix form, G = (A, B 7 ), 
where the payoffs are given as 


G = 




\(A n ^,Bl i) 


{Aln y , Bj' ) 


( 8 ) 


(-^-n x n.y j ^n x n. 


). 


For identical players the game is symmetric (A = B and n x = n y = n) and the game is well described by 
a single payoff matrix. An additional symmetry occurs in the partnership games iim (or ’’games with 


common interests 


271 ]) when the players share the utility equally, that is, if A = A 1 . 


For the two-player potential games we can introduce a potential matrix 

( V\\ ■ ■ ■ V\ n . 


V = 


(9) 


\V nx l ••• V r 


Tlx Tly _ 


where the matrix components can also be expressed via the use of a matrix product, namely, 


Vij — s X 2 • Vs yj , (19) 

that satisfies the conditions 

Vkj 1); — d /.! / d ij 1 

Vu - Vij = Bu — Bij , ( 11 ) 


where i, k = 1, ..., n x and j,l = 1, ..., n y . 

In order to visualize the relevance of the above constraints we introduce now the dynamical graph 
representation of games where each possible strategy profile is denoted by a node arranged in the same way 
as in the bi-matrix formalism [see Eq. ([5])]. 

Figure Q] illustrates the dynamical graph for a two-player 3x3 game where the edges of the graph connect 
those pairs of strategy profiles where only one of the players changes her strategy. This arrangement of nodes 
resembles a 3 x 3 square lattice with periodic boundary conditions. The payoff pairs and the elements of 
the potential matrix are indicated within the boxes representing strategy profiles. Notice, that the nodes of 
a column (or row) form a complete subgraph for any finite number of strategies. Along these edges we will 
consider the payoff variation of only that player who modifies her strategy. 

The existence of potential implies that the summarized potential variation (or payoff variation of the 
active player) between any two strategy profiles [e.g. along a given series of unilateral strategy changes 
from ( s X i,s V j ) to ( s x k,s y i )] is independent of the path connecting the initial and final strategy profiles. In 
other words, the sum of potential variation is zero along each loop of the dynamical graph. More precisely, 
according to Kirchhoff’s laws we should take into consideration only the independent loops, their number is 
the difference of the number of edges of the whole dynamical graph and those of its spanning tree [HU, |H^] . 

For the 3x3 games the dynamical graph has 9 nodes and 18 edges as plotted in Fig. |T] while the spanning 
tree has only 8 edges, e.g., those forming a shape of 3. The number of independent loops is 10 and the 
corresponding loops can be constructed by adding edges to the spanning tree consecutively. It is convenient 
to select the shortest loop containing the new edge. In the present case, however, the number of non-trivial 
independent loops is only 4 because along the three-edge loops (within a column or row) the conditions are 
satisfied as it is demonstrated by Eq. ©• 

Accordingly, in the present case there are only four independent (internal) four-edge loops that can be 
selected as denoted by dashed (red) circles in Fig. [TJ Kirchhoff’s law ensures that the conditions of the 
existence of potential are satisfied for all the possible loops if it is satisfied for each independent four-edge 
loop. 
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Figure 1: (Color online) Dynamical graph with payoff pairs and potential values for a two-player 3x3 game. Dashed (red) 
circles indicate a possible choice of the four independent four-edge loops to be considered for the existence of potential according 
to Kirchhoff’s law. 


For general treatment it is useful to discuss in detail the 2x2 sub-games where player x can use her 
strategy i or j (1 < i < j < n x ) while player y is restricted to select either her kth or Zth strategy 
(1 < k < l < n y ). The relevant payoff variations for this sub-game are illustrated in Fig. [2] A potential 



Figure 2: Payoff variations of the active player for a 2 x 2 sub-games along the directed four-edge loops of the dynamical graph. 
Within the boxes the upper parameters refer to the strategy profile and the corresponding payoffs are denoted in lower row. 


exists if the sums of the payoff variations of the active player along these four-edge loops become zero, that 
is, 

Bu — + Aji — An + — Bij + A^ ~ Ajk = 0 ( 12 ) 

for all possible 2x2 sub-games [27|. Conversely, all the sub-games of a potential game are potential games, 
too. 

For a two-player game with n x and n y strategies we can distinguish n x (n x — l)n y (n y — l)/4 two-strategy 
sub-games whereas the number of independent four-edge loops is significantly less, i.e. ( n x — 1 ){n y — 1), 
according to the Kirchhoff laws discussed above. The number of independent and relevant four-edge loops of 
the dynamical graph is reduced drastically for the symmetric matrix games. The reduction of the number of 
the relevant four-edge loops is related to the equivalence of payoff variations along the loops (i, k) (i, l) —> 
(j, 0 -t (j, k) ->• (i, k) and ( k,i ) (Z,i) ->• {l, j) (fc, j) (fc,i) if A = B. As a result, if condition m is 
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satisfied along a loop, then it is satisfied also for its counter-loop. 

Additionally, if A = B, then we can distinguish three types of four-edge loops. In the first case, both 
players can choose the same two strategies ( e.g ., i = k and j = l ) and then Eq. (fT^l) is always satisfied. 
That means that all the symmetric 2x2 games (and sub-games) are potential games. Section HOI is devoted 
to discuss the consequences of this inherent feature. 

In the second case the players have a common strategy and their second strategies are distinct. For 
example, the upper right loop, indicated by the (red) dashed circle in Fig. |T| represents such a situation 
when Eq. in becomes 

A 12 — A 21 + A 23 — A 32 + A 31 — A 13 = 0 . (13) 

Due to the above-mentioned symmetries this is the only criterion for the existence of potential in the set of 
symmetric 3x3 matrix games. 

The third type of the four-edge loops can appear only for n > 4 because here the players use different 
strategy pairs. For example, along the loop (1, 3) —> (1,4) —> (2,4) (2, 3) —> (1, 3) in a four-strategy game 

Eq. m obeys the following form: 

A 13 — A 31 + A 24 — A 42 + A 32 — A 23 + A 41 — A 41 = 0 . (14) 

The structure of Eqs. (fl3l) and (THl) implies a relationship between the existence of potential and the absence 
of cyclic components of the basis games (see Sect. ED- 

In summary, for an n x n matrix game, Kirchhoff’s laws give (n — l ) 2 conditions to be satisfied by the 
payoff matrix elements (Ay and By) for the potential games. If A = B then the number of independent 
conditions is reduced to (n — l)(n — 2)/2 in agreement with the fact that each symmetric 2 x 2 matrix game 
is a potential game and for the symmetric 3x3 potential games the nine payoff parameter Ay must satisfy 
only one simple linear relation given by (1131) . 

2.3. Flow graphs 

The flow graph is the directed version of the dynamical graph where each node represents a strategy 
profile and the directed edges connect those strategy pairs which differ in only one of the player’s strategies. 
The arrows on these edges point towards the strategy profile resulting in higher income for the player 
modifying her strategy unilaterally. Notice that here the (bold) arrows differ from those we used in Fig. [2] 
This approach limits the analysis to those cases where the preferences of these unilateral changes are defined. 

Some of the properties of potential games are straightforward consequences of well-known results of the 
theory of graphs, especially of directed graphs [H3|. The flow graphs are simple directed graphs that have 
neither parallel edges nor self-loops. Additionally, the flow graph of a potential game is free of directed loops 
as mentioned before. 

Figure [3] represents a flow graph for a three-strategy game. Due to the intimate relationship between the 
dynamical graphs and flow graphs here we use the same arrangement of strategy profiles as in Fig. [I] We 
emphasize, however, that for larger number of strategies (n x ,n y > 3) the structure resembling the square 
lattice with periodic boundary conditions is preserved whereas the graph should be extended by additional 
directed edges connecting any pairs located within the same column or row. Remarkably, directed loops 
cannot occur within one row (and column) of the flow graph of a matrix game because of the conditions 
defined by Eq. ©■ Furthermore, within each of these directed subgraphs there exists only one strategy 
providing the best income for the given player. As a strict Nash equilibrium is an optimal choice for both 
players, therefore the latter fact explains why the number of strict (pure) Nash equilibria is limited by 
min (n x ,n y ). 

Figure [3] shows the flow graph of a three-strategy game having only one strict Nash equilibrium that 
can be found by starting the search for any strategy profile and allowing the players to choose unilaterally 
a better strategy from their own point of view. After some steps the walk in this directed graph ends in a 
Nash equilibrium represented by a node without outgoing edges. Evidently, a similar flow graph holds for 
games where the payoffs are modified by smaller components that are not capable of reversing the directions 
of arrows. 
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Figure 3: Flow graph for a three-strategy matrix game. Labels of vertices/boxes refer to the strategy profiles and the directed 
edges indicate the preferred strategy change for one of the players. 


In the textbook of game theory the reader can find different methods for the determination of Nash 
equilibria. One of the standard methods is based on the iterated elimination of the strictly dominated 
strategies. In the notation of Sect. 12.11 strategy s x of player x is strictly dominated by strategy s' x if 
u x (s x , S- x ) < Ux(s' x , S- x ) Vs-x- Such a dominance can be easily recognized in a flow graph. For example, 
Fig. [3] represents a game where strategy 3 is dominated by strategy 2 as the arrows connecting the nodes of 
the second and third columns are pointed to left. Notice, that for this special case strategy 1 dominates the 
other strategies for both players. 

As rational players do not play dominated strategies, these strategies can be eliminated. Once we 
have eliminated a dominated strategy it can happen that another strategy for one of the players becomes 
a dominated strategy. For the games illustrated in Fig. [3] only the (1,1) strategy profile survives after 
repeating this procedure. 

The above method will not eliminate rows and columns including a strict Nash equilibrium. The iterated 
elimination of dominated strategies simplifies the game and shrinks the flow graph. In many cases the 
simplified model may be equivalent to a well-known game ( e.g ., coordination game). 

The number of Nash equilibria is a fundamental question in the theory of games. In the literature of 
graph theory 55, HiJ the reader can find the discussion of finite directed acyclic graphs, i.e. finite directed 
graphs that have no directed loops. Such graphs always have at least one node without outgoing edges, 
otherwise a directed walk of arbitrary length could be constructed in the finite graph, as we could always 
leave any vertices of the loop, and by the finiteness of the graph the walk would contain a cycle, contradicting 
the assumption of acyclicity. 

The latter statement specifies the Nash theorem for the case of potential games. Namely, at least one 
strict pure Nash equilibrium exists for the potential games. By reversing the direction of the edges, a similar 
argument shows, that a finite acyclic directed graph also has at least one node without incoming edges. The 
corresponding node is represented by the strategy profile (3, 3) in Fig. [3] 

Remarkably, the determination of flow graph for any matrix game can help us find the Nash equilibria as 
we only need to identify the strategy profiles without outgoing edge(s). In the rest of this work we discuss 
more complex flow graphs corresponding to multi-player games where many Nash equilibria can exist. 


2-4- General features of the potential 

The linear relationship between the individual payoff variation and potential variation implies general 
features. First we emphasize that a game as well as its potential is not influenced if all payoff elements 
(including Ay, By, and Vy-) are increased with a constant value because the player decisions depend only 
on the payoff differences. 
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If a bi-matrix game G = (A, B T ) has a potential V then the game G' = (aA, aB 2 ) also has a potential 
given as V' = aV. If a > 0 then the multiplication of the payoff elements can be interpreted as a choice 
of new unit. For negative a, however, the game and its Nash equilibria can change drastically, whereas the 
existence of potential remains valid. 

In addition, if we have two potential games, G = (A, B T ) and G' = (A', B ,T ), with potentials V and V' 
then the game obtained by the linear combination of payoffs G * = qG + f3G' = {aA + 0A', aB T + /3B ,T ) 
is also a potential game with V* = aV + /3V'. 

One can easily check that for symmetric matrix games (A = B) the potential is a symmetric matrix, 
i.e., V = V T . This feature is the direct consequence of the fact that players exchange their payoff if 
they exchange their strategies. An additional symmetry occurs for the partnership games (A = B and 
A = A 2 ) when V = A. For such situations the players share their income equally, resembling the fraternal 
or egalitarian behavior [US, [HU. In that case, the individual and common interests coincide and many 
intriguing events (e.g., social dilemmas) are dropped. 

Within matrix games we can distinguish games with self- and cross-dependent payoffs j^j] . For the cross¬ 
dependent payoffs the income of each player depends only on the co-player’s strategy, that is the columns 
of the payoff matrices are composed of uniform values ( 7 j and 5i for j = 1 ,... ,n y and i = 1 ,..., n x ) as: 



( 7i • 

• 7 n y \ 


(S 1 ■ 

&n x \ 

A = A< cr > = 

\7i • 

■ ln v / 

, B = B < ' cr ) = 

U • 

firix J 


In these cases the unilateral strategy changes are not motivated by receiving a higher individual payoff 
and the corresponding potential is constant (i.e., we can choose V( cr ' = 0). On the contrary, for the 
self-dependent payoffs the rows of the payoff matrices are constant, namely, 



l 7i ’ 

■ 7i \ ( h ■ 

■ «1\ 


A = A (s) = 


: , B = B (s) = I : 

! ’ 

(16) 


\"fn x ■ 

7 n x J V S ny 

• SnJ 



and the components of the potential matrix obey the following form: 

^S S) =7 i + Sj, (17) 

or in matrix notation 

V (s) = A + B t . (18) 

The linear relationship between the potential and the payoffs offers the possible use of decomposition 
when the matrix games are composed of elementary games reflecting basic symmetries as detailed in the 
next section. 


3. DECOMPOSITION OF TWO-PLAYER GAMES 


In the literature of game theory the concept of decomposition is used in different ways. For most of the 
cases the games with many payoff parameters are built from simpler games characterized by a significantly 
less number of parameters [see the works by Szep and Forgo [56 


Kleinberg and Weiss 57| . Sandholm 


[53 . and Candogan et al. [Q]. This gives us a deeper insight into the general properties of interactions 
described by the tools of game theory. Relationship between payoffiSj potential, flow graphs, and other types 
of interactions are discussed in the papers by Candogan et al. [ 2 ^, [58|, Hwang and Rey-Bellet [30] who used 


a different terminology. Now we follow the concept introduced recently in the papers by Szabo et al. 31 


which is based on the introduction of orthogonal elementary basis games representing proper features, e.g., 
games with self- or cross-dependent payoffs. Our analysis will be focused on those decompositions that help 
us identify the components prohibiting the existence of potential. This way of decomposition is consistent 
with the stability analysis based on the systematic investigations of two-strategy sub-games [6(JEl[. 
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3.1. Decomposition of symmetric matrix games 

Symmetric matrix games are defined by the Aj :j elements of the payoff matrices (B = A and n x = n y = 
n). In the spirit of decomposition, the payoff matrix is a linear combination of four matrices: 


A = 


An A12 

A 21 A 22 


= A 1 


1 0 
0 0 


A 1 


0 1 
0 0 


Ao 


0 0 
1 0 


Ao 


0 0 
0 1 


(19) 


with the suitable coefficients. For a linear arrangement of the matrix components the above four matrices 
form a set of orthogonal ” basis vectors” that are called basis matrices or basis games henceforth. We can 
choose, however, another set of orthogonal basis matrices representing fundamentally different games that 
involve the relevant symmetries. For example, a new set of orthogonal basis matrices can be defined as 


f(l) = 


1 1 
1 1 


, f(2) = 


1 -1 

1 -1 


f (3) = 


1 1 
-1 -1 


, f(4) = 


1 -1 

-1 1 


( 20 ) 


These basis games were used in many examples studied previously [see e.g. 
defined by introducing the scalar product of two n x n matrices as 


H S3]- The orthogonality is 


A ■ B — B • A — 'y ' Aij Bij , 
i,j =1 


( 21 ) 


that is zero if A and B are orthogonal to each other. Notice that the basis games given by 
CUD form two sets of orthogonal basis matrices as f(m) • f(m') = 0 if m ^ m' ( m,m' = 
otherwise f(m) • f(m) = M{m) > 0. 

Using the orthogonal basis matrices one can express any n x n payoff matrix as 

n 2 

A = ^2, ct(m)f(m) , (22) 

m= 1 

where the coefficients a (to) are given by the expression 

1 1 71 

a(w) = %) A ' f(m) = A/fm) (23) 

The first basis matrix f (1) given by (CiUl) defines a game with constant payoffs. In fact this is the irrelevant 
component of the game that does not influence the decision of selfish players and a(l) [defined by Eq. (l22l) ] 
is equal to the average value of A^. The linear combinations of the first and second components {e.g., 
A = A^ cr ) = a(l)f(l) + a(2)f(2)) describe all games with cross-dependent payoffs. Similarly, the matrix 
A = A( s ) = a(l)f(l) + a(3)f(3) corresponds to a game with self-dependent payoffs. 

In game theory f(4) represents the coordination game when players are enforced to choose the same 
strategy. At the same time, the game with payoff matrix A = — f(4) is an anti-coordination game where 
the best results can be achieved by the players if they choose opposite strategies. Thus the sign of a(4) 
defines whether coordination or anti-coordination type interaction is built into the game and its strength is 
measured by |a(4)|. 

In the light of the above results, a symmetric two-strategy game is composed of three types of components 
representing the self-dependent and cross-dependent payoffs, and the coordination type interaction. Due to 
the general features discussed in Sec. 12.41 only two of these components contribute to the potential matrix 
that obeys the following form: 

V = a(3)[f(3) + f T (3)]+a(4)f(4) (24) 

omitting the arbitrary constant [proportional to f(l)]. 

The mentioned general features are preserved when the number of strategies is increased, whereas a 
fundamentally new type of interaction emerges if n > 2 (or even if B 7 ^ A). This type of interaction 
represents the cyclic dominance that prevents the existence of potential. 


Eqs. (mill and 
1 ,..., n 2 ) and 
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For the illustration what happens for n > 2, we briefly discuss now a possible decomposition of the 
symmetric 3x3 matrix games. A more detailed analysis is available in a recent paper by Szabo et al. [31] 
who suggested the introduction of a set basis matrices representing the two-dimensional Fourier components 
of A. Instead of it, now we use another notation based on the dyadic products jHtjj]. In this notation the 
basis matrices g(m) (m = 1,..., 9) are expressed with the help of dyadic products e(k) 0 e(l) ( k , l = 1,2,3) 
of the following three three-dimensional orthogonal basis vectors 


e (l) = 



e(3) = | -1 


(25) 


The relationship between the matrix label m and the vector labels k and l is given by the following definitions: 

(26) 


g(l) = 

g(2) = 

g(3) = 

g( 4 ) = 
g(5) = 

g(6) = 

g(7) = 

g(8) = 

g(9) = 



- [e(2) 0 e(3) + e(3) <g> e(2)] = 


- [e(2) 0 e(3) - e(3) 0 e(2)] = 



(27) 

(28) 

(29) 

(30) 

(31) 

(32) 

(33) 

(34) 


Accordingly, the 3x3 payoff matrix is expressed by the linear combinations of these g(m)s as 

9 

A = > 

m —1 

where the coefficients j3(m) are given by the expression 

Q 


^ (m) = AfR A ' g(w) = ATR £ 


(35) 


(36) 


j,k= 1 
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and the normalization factors are defined as above, that is, A f(m) = g(m) • g(m) . 

The present set of basis vectors reflects similar features and symmetries we found for the 2x2 payoff 
matrices. Namely, the first component g(l) defines a constant (or average) contribution to the payoffs and 
for later convenience its contribution will be denoted for n > 2 as 

A (av > = /3(l)g(l) (37) 

where g(l) denotes the all-ones matrix and 

1 n 

m = - 2 £ ( 38 ) 

i,3 =1 


The linear combination of the first three terms [A^ cr ) = /3(l)g(l) + /3(2)g(2) +/3(3)g(3)] spans the whole 
set of symmetric 3x3 games with cross-dependent payoffs. On the other hand, the self-dependent payoff 
component of a game A can be given as A (s ) = /3(l)g(l) + /3(4)g(4) + /3(5)g(5) where values of /3( 1), /3( 4), 
and /3(5) are given by Eqs. (l36l) . Straightforward calculations yield that these portions of payoffs obey a 
simple form: 


A« = 


'71 

7i 

71 \ 

72 

72 

72 

^73 

73 

73/ 

in the ith 

(i = 


(39) 



3 


£ Aij 


(40) 


The values of 7 are equivalent to the average income of player x when she chooses her ith strategy, while 
her co-player follows a mixed strategy by choosing the three options with equal probabilities. For games 
with more than one Nash equilibria the concept of risk dominance, suggested by Harsanyi and Selten [40]. 
gives us an additional criterion to select one of them. The risk dominance dictates to choose those strategies 
that provide the highest 7 , value. Accordingly, A^ s ) quantifies the risk dominance, too. Figure [3] illustrates 
the flow graph of A( s ) for 71 > 72 > 73 - For games A^ a potential exists with a potential matrix 
= A^ + A( s > t as it is detailed in Sect. 12.41 

Evidently, one can perform a similar calculation for the evaluation of the component A^ cr ^ = /3(l)g(l) + 
/3(2)g(2) +/3(3)g(3) with cross-dependent payoffs that can be expressed as A^'p = Aij /3. These terms 

do not contribute to the potential matrix (V( cr ) = 0). 

The component g(7) [see (1421) ] corresponds to an elementary coordination type game where the players 
achieve the best if both choose either the first or the second strategy simultaneously. Due to its relevance 
this payoff matrix is denoted henceforth as d(l,2). Additionally, we can introduce further elementary 
coordination games d(p, q) defined as 


dij(p,q) 


1 , if i = j=p, 

1 , if i = j=q, 

—1, if i = p and j = q , 
— 1 , if i = q and j = p , 
0 , otherwise, 


(41) 


even for n > 3 when 1 < p < q < n that describes similar relationship between the strategies p and q (for 
simplicity the n-dependence is not denoted). Notice that d(l, 3) can be obtained from d(l, 2) by exchanging 
the second and third rows and columns simultaneously. Although the matrices d(l,2), d(l,3), and d(2,3) 
are not orthogonal to each other, these basis games span the three-dimensional subspaces of the coordination 
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type interactions defined as A^ coord ) = /3(6)g(6) + /3(7)g(7) + /3(8)g(8). Applying the expressions © and 
(IMl) one can easily check that, for example, 


d(l, 3) 



\ [g(6) + g(7) - 2g(8)]. 


(42) 


The anti-coordination type interactions between any strategy pairs (p < q ) are also located along the 
half-lines (/3d(p, q) with /3 < 0) in the space of A( coord ). For example, the game with a payoff matrix 
A(c°° rd ) = — g(7) corresponds to a situation where the players have two equivalent Nash equilibria, namely 
the strategy pairs (1, 2) and (2,1). 

The coordination type games include curious cases when the game has three equivalent Nash equilibria. 
For instance, the strategy pairs (1,1), (2, 2), and (3, 3) are equivalent Nash equilibria in a game defined as 
A = d(l,2) + d(l,3) + d(2,3). If one exchanges the first and second rows and columns simultaneously in 
the latter payoff matrix, then the resultant game is also a coordination type game with three equivalent 
Nash equilibria: (1,2), (2,1), and (3,3). 

Notice that A( coord ) T = A^ coord ) indicating the coincidence of individual and common interests within 
this type of games and the existence of a potential matrix that is identical to the symmetric payoff matrix 


[■y(coord) _ ^(coord)j 

The last component (g(9) = A (rps )) represents the well investigated rock-paper-scissors game that has 
only one mixed Nash equilibrium when the three strategies are selected with the same probability. The 
rock-paper-scissors games exemplify those systems from ecology to non-equilibrium physics, where three 


states (strategies or species) cyclically dominate each other [62], |4|, |63[ [6J, .6 


The rock-paper-scissors basis game (g(9)) has no potential because of the existence of directed loops in 
the flow graph shown in Fig. [4] For the illustration of the three-fold symmetries within this flow graph, here 
the strategy profiles (nodes) are rearranged. This flow graph illustrates clearly the presence of a directed 
six-edge loop along the hexagonal periphery where the thick directed edges represent the ” best response”. 
Besides it, this flow graph contains six directed four-edge loops [ see e.g ., the loop (1, 2) (1,3) (2, 3) —► 



Figure 4: Flow graph of the rock-scissors-paper game. Thick edges represent the best responses. 

(2,2) -A (1,2)]. 

As potential exists for all the linear combinations of the first eight basis games, therefore the absence of 
the rock-paper-scissors component (/3(9) = 0 ) can be interpreted as a necessary condition for the existence 
of potential. The latter condition, more precisely the mathematical expression of A-A^ rps ^ = 0, is equivalent 
to m we derived previously with the application of the Kirchhoff laws. 
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The payoff matrix g(9) = A*f' ps ) can be considered as the adjacency matrix of a simple directed graph 
(shown in Fig. [5]) that illustrates graphically the cyclic dominance among the three strategies represented 
by the nodes of this dominance graph. In order to avoid confusion with the dynamical and flow graphs here 



Figure 5: The adjacency matrix of this three-node directed graph is equivalent to payoff matrix of the rock-paper-scissors game. 
The directed three-edge loop refers to cyclic dominance among the three strategies denoted by labeled black bullets. 


we used a third way for the graphical illustration of a simple directed graph. The conventions of the arrow 
directions in the flow and dominance graphs become identical if the adjacency matrix elements of a directed 
graph are defined as Ay = —Aji = 1 if there is a directed link from the strategy (node) j to i and otherwise 
Aij = 0. In words, the dominance graph in Fig. [5] corresponds to a game where strategy 1 dominates 2, 2 
dominates 3, and 3 dominates 1. This direction convention ensures that the arrows will show the direction 
of probability current in evolutionary systems discussed later. 

It is noteworthy, that the cyclic dominance graph is analogous to the cyclic food web used in ecological 
models to describe the predation-prey relationship. The adjacency matrix of the tree-like dominance graph 
(trophic food web) can be related to some anti-symmetric components in the sub-space of games with the 
self- and cross-dependent payoffs. 

Evidently, one can choose another set of basis matrices with features suitable for the problems to be 
studied. For example we can separate the symmetric and anti-symmetric parts of the self- and cross¬ 
dependent payoffs which allow the analysis of the zero-sum components or the strength of social dilemmas. 

Most of the general properties of the symmetric three-strategy games are preserved for n > 3. The 
decomposition of a symmetric matrix game into elementary components becomes impressive for n = 2 k (k is 


an integer) when the columns of the Walsh-Hadamard matrices [66j are used as basis vectors in the dyadic 
decomposition Q. For this choice both the basis vectors and the derived basis matrices are composed 
of only +ls and —Is that simplify the calculations [for an example see Eq. (12011 ]. These calculations have 
confirmed the relevance of the mentioned four types of interactions for n = 4 [H^]. Similar results are 
concluded by Candogan et al. [29|, Hwang and Rey-Bellet [30] who studied the decomposition of ?r-strategy 
games without introducing a definite set of basis games. 

According to the analyses mentioned, the games with self-dependent payoffs are spanned by n dyadic 
products, as e(l) g e(j) for j = 1,..., n , of an orthogonal set of basis vectors (e(l),..., e(n)) when the 
first one is composed of only +ls, that is, e,;(l) = 1 for i = 1,..., n. In this notation the basis games with 
cross-dependent payoffs are given as e(j) g e(l). The all-ones matrix (g(l) = e(l) g e(l)) belongs to both 
types and plays the role of irrelevant term. Its coefficient /3(1) quantifies the average payoff. The direct 
pair interactions are missing within their unified parameter space spanned by the (2 n — 1) orthogonal basis 
matrices that may even be divided into the sum of symmetric (e(l) g e(j) + e(j) g e(l)) and anti-symmetric 
(e(l) g e(j) — e(j) g e(l)) basis matrices. Thus the unified parameter space of A (cr ) and A*®) can also be 
spanned by the mentioned n symmetric and (n — 1) anti-symmetric basis games. 

The rest of the symmetric parameter space of games is equivalent to the linear combinations of the 
additional symmetric basis matrices defined by the dyadic products as e(fc) g> e(j) + e(j) g) e(fc) if 1 < k < 
j < n. Hwang and Rey-Bellet [30] have shown that this set of games (denoted as A (coord )) are also spanned 
by the n(n — l)/2 d(p, q) matrices defined by Eq. (Hill (for 1 < p < q < n). Notice that all the mentioned 
basis matrices of A (coord ) are orthogonal to both A^ and A (cr ) because the sums of the matrix elements 
are zero in the rows and columns separately. The number of these orthogonal basis games, as well as the 
number of the independent d (p,q) matrices, is n{n — 1 )/2. 

The rest of the anti-symmetric parameter space of the payoff matrix is spanned by the orthogonal basis 
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matrices e(k) <8> e(j) — e(j) <g) e(fc) where 1 < k < j < n. The number of these orthogonal basis games is 
(n — l)(n — 2)/2 and their linear combinations are also denoted by because of the presence of cyclic 

dominance. 

Calculations [59| for n = 4 have indicated that A( cycl ) is spanned by three orthogonal basis games 
defined by the adjacency matrices of the four-edge directed graphs plotted in Fig. [6] Additionally, A( cycl ) 



Figure 6: The adjacency matrices of these four-edge directed graphs are orthogonal to each other and represent the three cyclic 
basis games for n = 4. 

contains games where the cyclic dominance is limited to three strategies as in the rock-paper-scissors game. 
For example, the sum of these three cyclic orthogonal basis games is equivalent to a rock-paper-scissors 
type sub-game (with strategies 2, 3, and 4). The knowledge of these three cyclic basis games simplifies the 
identification of potential games, because it is enough to check whether the scalar products of A by the 
three cyclic basis games become zero simultaneously, or not. Potential exists if all these scalar products are 
zero. 

For n > 4, however, determining the existence of potential becomes difficult because of the large number 
of the possible four-edge and three-edge loops illustrated in Fig. [7] It is emphasized, that the orthogonality 

(a) O W O 

© 0 O 0 



Figure 7: Directed graphs with n nodes and a single three- (a) or four-edge (b) directed loop. 

between A and the adjacency matrices of these types of the directed graphs reproduces the same conditions 
for the existence of potential that we have derived in Sec. 12.21 The application of the Kirchhoff laws gives 
us a method to select (n — l)(n — 2)/2 independent loops to be checked. 

In sum, the payoff matrix of a two-player symmetric game can be decomposed into the sum of four 
classes of payoff matrices representing the self-dependent, the cross-dependent, the coordination, and the 
cyclic dominance type interactions. More precisely, 

A = A< 8 > + A( cr) - A< av) + A( coord ) + A (cycl ), (43) 

that takes into consideration that A ( av ) is involved in A^ s ) and A( cr ). Potential exists if A( cycl ) = 0 and the 
potential matrix can be given as 

V = A^ + A (s)t + A^ coord \ (44) 


3.2. Decomposition of two-strategy bi-matrix games 

The decomposition of the bi-matrix games requires the straightforward extension of the scalar product 
of two bi-matrix games analogously to the traditional scalar product of two vectors. As the bi-matrix game 
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G is described by the elements of two matrices (A and B), therefore the scalar product of the games G and 
G' (determined by (A' and B')) is defined as 


G G'=^ [AijA' i:j + BijB'ij]. 

i,j =1 


(45) 


The decomposition of the two-strategy games becomes impressive if we use the following orthogonal basis 
games: 


f'(i) 

f'(2) 

f'(3) 

f'(4) 

f'(5) 

f'(6) 

f'(7) 

f'(8) 


(( 1, 1) ( 1, 1)\ 
V( 1, 1) ( 1, 1 )) 
(( 1, 1) (-1, 1)\ 
l( i.-i) 

(( i, i) ( i,-iA 

V(-i. i) (- 1 ,- 1 ); 

(( i, i) (-i,-i)A 
V(-i.-i) ( i, i); 
({ i,-i) ( i,-d\ 

v( 1 ,- 1 ) ( i,-i); 

/( i,-i) (-i,-D\ 

V( i, i) (-i, i); 
(( i,-i) ( i, in 
V(-i.-i) (-i, i); 

(( i,-i) (-i, i)\ 

v(—1, 1) ( 1,-1)) 


(46) 

(47) 

(48) 

(49) 

(50) 

(51) 

(52) 

(53) 


that are composed of +ls and —Is. Using these orthogonal basis games, all the two-strategy bi-matrix 
games can be given as 

8 

G = a'{m)i'{rn) (54) 

m =1 

where the coefficients a'(n ) are expressed as before, i.e., 

a'(m) = -^G • f'(m). (55) 

8 


Notice that the first four basis games span the parameter space of the symmetric games when B = A, 
that have been discussed previously. In other words, f'(l),..., f'(4) denote the bi-matrix version of the 
symmetric two-strategy orthogonal basis games given by Eqs. ®. 

The linear combinations of the additional second four basis games describe the anti-synrmetric games 
where B = —A. The resulting four basis games [U(5),..., f'(8)] are obtained from the first four ones by 
reversing the sign of payoffs received by the second player. All these basis games represent some properties. 
For example, the coefficient of U(5) quantifies the difference in the average payoffs between the players x 
and y. 

The games with self-dependent payoffs are given as 

G (s) = a'(l)f'(l) + c/(3)f'(3) + c/(5)f'(5) + c/(7)f'(7), (56) 

and the resultant payoff matrices A^ s ^ and B^ s ' are composed of uniform rows with payoff parameters 
representing averages values as for the symmetric games ( e.gA[^ = (Aji+Ai 2 )/2 and = {Bn+Ba)/2). 
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Similar expressions define the components with cross-dependent payoffs: 

= a'(l)f'(l) + a'(2)f'(2) + e*'(5)f'(5) + a'(6)f'(6). (57) 

The basis game f'(4) defines the symmetric coordination type interaction with payoffs A (coord ) = 
g (coord) _ a '^4)f(4) where f(4) is defined by Eq. (15U1) . 

The bi-matrix game f'(8) is equivalent to the matching pennies game, which is a well-known zero-sum 
game where the first player wins a payoff unit from the co-player if they both choose the same strategy. For 
opposite choices the second player wins 1 from the first player. This game has a mixed Nash equilibrium, 
where the players choose their strategy at random. 

The flow graph of the matching pennies game (see Fig. [5J illustrates the appearance of a directed loop 
because for each strategy profile the unsatisfied player can increase her own payoff by 2 by reversing her 
strategy. As a result, if in an iterated game the randomly selected player is allowed to modify her own 



Figure 8: For the matching pennies game one of the players always wishes to reverse her strategy in a way maintaining cyclic 
strategy modifications. 


strategy unilaterally then the actual strategy profile will evolve cyclically. Such situations were reported 
previously by van Valen 67, 68] who suggested the concept of Red Queen mechanism to explain biological 


evolution via a constant arm race between co-evolving species. The interactions between buyers and sellers 
[69j |. property owners and criminals 0, or conformists and rebels |Yl| exhibit similar features. For very 
recent investigations of the Red Queen mechanism at the level of population dynamics we suggest reading 
the papers by Sardanyes and Sole (72| and Juul et al. [73]. The cyclic variation in the strategy profiles has 
been observed by Xu et al. jY4j in human experiments. 

Evidently, potential does not exist for the matching pennies game and also for those two-strategy bi¬ 
matrix games that include this component. It is needless to emphasize that the bi-matrix game G = G*- pot ) 
is a potential game, if a'(8) = G • g(8) = 0, and the latter existence criterion for the potential is equivalent 
mathematically to Eq. CE2D for the particular case when i = k = 1 and j = l = 2. 

For the rest of (seven) basis games the potential of a non-symmetric bi-matrix game can be easily 
determined. According to this calculation, the potential matrix V( bmg ) can be built up only from three 
components as 

V (bmg) = c/(3)[f(3) + f T (3)] + a'(4)f(4) + a'(7)[f(3) - f T (3)]. (58) 


It is worth mentioning that the game G = G( pot ) +a'(8)f / (8) can even be an ordinal potential game with 
the potential of G( pot ) [given by Eqs. (1551) ] if |a'(8)| does not exceed a threshold value dependent on the 
payoff differences in G*- 150 ^, more precisely, if the contributions of the matching pennies game cannot reverse 
the arrow directions in the flow graph that ensures the existence of at least one pure Nash equilibrium. 
Conversely, if the cyclic component dominates the system behavior, then the game has only a mixed Nash 
equilibrium until G^ pot ^ is weak enough to reverse the flow direction along at least one edge in the flow 
graph. 


3.3. Properties of two-player two-strategy games 

In the previous sections we have shown a general method to evaluate the potential matrix in the absence 
of cyclic (or matching pennies) components. This method is based on the concept of decomposition and 
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exploited the symmetries simplifying the calculations. Now we introduce another method and notation to 
discuss the general properties of the two-player two-strategy games. 

For a deeper discussion of the symmetric two-strategy games (A = B) we use the traditional notation 
introduced for the investigation of social dilemmas [0, 0,0 where the strategies are called defection and 
cooperation (in short s x i, s y \ = D and s x 2 , s y 2 = C). Within this terminology both players receive P (pun¬ 
ishment) or R (reward) for mutual defection or cooperation meanwhile the defector receives T (temptation) 
and her cooperative co-player gets S (sucker’s payoff) if they choose opposite strategies. For the case of 
prisoner’s dilemma (S < P < R < T) both selfish players are enforced to choose defection (this is the pure 
Nash equilibrium) because the cooperation is a dominated strategy, that is, the players receive higher payoff 
when choosing D irrespectively of the co-players decision. The curiosity of the prisoner’s dilemma game is 
that the players’ rational (selfish) choices yield a suboptimal outcome. 

There are two other (weaker) social dilemmas, namely the hawk-dove [ 3 } (called also snowdrift or 
chicken [78j) game (with a payoff rank of P < S < R <T) and the stag hunt game (S < P < T < R). Most 
of the previous analyses are performed when setting P = 0 and R = 1 without loss of generality. 

In the terminology of social dilemmas the bi-matrix of the symmetric two-strategy game is usually given 


as 


G = G( sd ) 


(0,0) (T,S)\ 

(S,T) (1,1) J- 


(59) 


For these payoff parameters the three social dilemmas are positioned on the S — T plane as illustrated in 
Fig. [H] where four regions are distinguished. 



Figure 9: (Color online) Four types of symmetric two-strategy games can be distinguished as a function of T and S for the 
notation of social dilemmas when the black and white bullets represent D and C strategies. The four segments are characteristic 
to games abbreviated as: H = harmony; HD = hawk-dove; SH = stag hunt; and PD = prisoner’s dilemma. The pure Nash 
equilibria are represented by nodes without outgoing edges in the flow graphs. 


Within each region the two-person game can be characterized by the same flow graph denoted. It is 
emphasized that these notations of strategies and payoffs are related to the case of the prisoner’s dilemma 
game. For other types of games one can find more suitable names for the strategies and payoffs, e.g ., for the 
hawk-dove game the players are to choose between the aggressive and peaceful (conflict avoiding) behaviors 
and the resultant payoff can be well described by using the terms of equal share of benefit, exploitation, and 
cost of injury. Anyway, the hawk-dove game has an additional mixed Nash equilibrium when the players 
choose one of their two strategies with a probability dependent on the payoffs. 


19 
































The absence of additional types of flow graphs is related to the symmetries that prescribe equivalent 
payoff variations for the active player deviating from the ( D,D ) or (C,C) strategy profiles. Due to these 
symmetries all the symmetric two-strategy games are potential games and the potential matrix is formally 
given by Eq. m with the decomposition of the payoff matrix 

A = A( sd ) = (° ^ (60) 


into the sum of four components [as defined by Eas. E^l) and m) with coefficients 


*(!) = 


T + S + l 


*( 2 ) = 


S —T — 1 


*(3) = 


T-S-l 


K4) = 


1 — T — S' 


According to Eq. (HH1) the potential obeys the following form: 


y( sd ) 





(61) 


(62) 


This expression will be used in Sec. 16. ll to illuminate the relationship between the Ising models and the multi¬ 
agent evolutionary games, if the interaction is defined by symmetric two-strategy games. It is emphasized, 
however, that the above potential matrix can be expressed by an even simpler expression, 


y(sd) 


0 S \ 

S 1 -T + SJ ’ 


(63) 


that is obtained from (El by adding a suitable constant to Vij s. The second version of the potential can 
be evaluated by using a simple algorithm. First we choose V\\ = An, subsequently V \2 = V 21 = A 21 , and 
finally V 22 = V 12 + A 2 2 — A 12 in the spirit of Eqs. ©■ This method can be adopted for larger number of 
strategies if the existence of potential is already justified. The resultant potential matrices are convenient 
when deriving phase diagrams for some types of multi-agent evolutionary games. 

It is worth mentioning that the values of the potential are equal for the two pure Nash equilibria [(C, D) 
and ( D,C)\ within the region of the hawk-dove game (see Fig. [9]). On the contrary, the potential values 
differ between the (D,D) and (C, C) Nash equilibria within the region of stag hunt games. To be more 
quantitative, (C, C ) is the preferred Nash equilibrium if S > (T — 1) when T < 1 and the strategy pair ( D, D) 
is preferred if S < min [(T — 1), 0]. The location of the preferred Nash equilibria is illustrated graphically in 
the next section (see Fig. [TO where it is contrasted with the case of collaborating players. 

Within the region of the hawk-dove (or anti-coordination) games the equivalence between (C, D) and 
( D,C ) Nash equilibria is evidently broken if A ^ B. For these bi-matrix games the potential ([58]) obeys 
the form: 

V <lm!, =°'(4)(_; ~[)+«'<3>(o l), (64) 

that we use later when drawing a parallel between the multi-agent evolutionary games and Ising models. 
Now we emphasize that the first components in Eqs. E2D and (1M1) measure the strength of coordination 
type interaction, the second ones quantify the potential difference between the (1,1) and (2, 2) strategy pairs 
whereas the third term distinguishes the potential value between the states (1, 2) and (2,1). 


3-4- Fraternal collaboration versus individualism 

In traditional game theory [lj the formation of coalition is permitted for the so-called cooperative games. 
Some of the essential elements of the games are dropped if the players are allowed to agree in how they wish 
to increase the sum of their income for two-player potential games. On the other hand, such situations can 
occur frequently in our everyday life. Thus we analyze now what happens if fraternal players try to maximize 
their total income shared equally. The results can serve as a reference when comparing actions motivated by 
individual or common interest. For this goal we introduce a reevaluated payoff matrix, A! fr ) = (A + A r )/2, 
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defining the individual payoffs for symmetric games with n strategies. For these games the individual and 
common interests coincide and the preferred Nash equilibrium may differ from those proposed for a potential 
game A. 

The equal share of payoffs eliminates the anti-symmetric components of A and the corresponding po¬ 
tential matrix is given as V( fr ) = A^ fr ). In the two-strategy social dilemma notation discussed above 


V( fr ) = 


0 (T + S)/2\ 

(T + S )/2 I ) 


(65) 


which predicts two distinct regions on the T — S plane when considering the preferred Nash equilibrium 
that may be either ( C,C ) (if (T + S) < 2) or the (C,D) or (D,C) strategy profiles (if (T + S) > 2). The 
latter results differ strikingly from those we obtained above as it is illustrated in Fig. [TUI The differences in 
the location of the maximum matrix elements of V [see Eq. ([63]) ] and V (fr ) imply social dilemmas for the 
two-strategy games. 



Figure 10: Preferred Nash equilibria for selfish (left) and fraternal (right) players. Dash-dotted lines divide the S — T plane 
into four regions nominated in Fig. [9] 

A similar comparison can be performed for other potential games even for n > 2. In that case we use 
the notation of Eq. C51) when potential exists if A^ cyc1 ^ = 0 and the resulting potential matrix is determined 
by only two components (A* s ^ and A* coord ^) as expressed in (HU). For these games the effective payoff and 
potential matrices for fraternal players involve the contribution of A^ cr \ too, as 

y( fr ) = A( fr ) = A^ coord - ) — A^ av ) + -[A^ + A^ T -(- A^ cr ^ + A^ cr ^ r ]. (66) 

We face ’’social dilemma” if the preferred Nash equilibria are distinct, that is, when (i,j) where 

max [14;] = Vij and max [ V fc ; '] = V-}-,’. The conflicting situations involve cases when the selfishness results 
in lower payoffs for both players (as it happens for the prisoner’s dilemma) or also those ones when the loss 
of the sucker exceeds the extra profit of the winner. The discrepancy between the individual and common 
interests can be quantified by 


V - V (fr) = i[A (s) + A (s)t - A (cr > - A (cr)T ] (67) 

depending only on the self- and cross-dependent components. More precisely, V — V^ fr ^ is the potential of 
the antisymmetric portion of the self- and cross-dependent payoff components defined as (A (s ) + Ai cr ^ — 
\( S ) T — A( cr ) T )/2 = A — A 7 . Accordingly, no dilemmas arise if A ^ = A^ cr T. In the light of this result 
for a symmetric n-strategy potential game we have (n — 1) independent parameters to control the absence 
or presence of the social dilemma in the above sense. 

The social dilemmas and their consequences are preserved for spatial multi-agent evolutionary games with 
equivalent players if logit rules control the strategy updates as it is demonstrated by Szabo and Szolnoki 
[79( 1 who studied cases of non-equal sharing for n = 2. 
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4. MULTI-PLAYER POTENTIAL GAMES 


Up to this point we have discussed only two-player games. Evidently the sum or any linear combinations 
of independent two-player potential games are also potential games as the original conditions © are satisfied. 
This situation remains valid even for the cases when some of the players participate in more than one potential 
games using the same strategies in each action. Remarkably, if a player uses two (or more) strategies against 
her co-players then she can be considered as two (or more) players. In many biological and economical 
systems the strategy represents the type of species or behavior that is the same for all the interactions for 
a given individual. For human systems, however, one can assume that the players use different strategies 
in the games they participate. Furthermore, human players can use mixed (stochastic) strategies in the 
repeated games. 

Henceforth we assume that the players use pure strategies expressed by unit vectors [see Eq. ([6]) ] for 
each player x = 1 ,.N where N denotes the number of players. In the present game theoretical models 
the players wish to maximize their own payoff by consecutive strategy refreshments disregarding the payoff 
variation for all the other co-players. If the players play two-person potential games with a subset of the 
other players, then their individual motivation can be well quantified by the variation of potential 

U(s) = ^2 s x -V( x ,y)s y (68) 

{x,y) 


where the summation runs over all interacting pairs (denoted as ( x,y )) defined by the edges of the connec¬ 
tivity graph where each node represents a player. 

In real social systems the pair potential may be different for all interacting pairs as indicated by the 
expression V(:r, y). In most of the cases, however, we assume that the players are equivalent and the games 
are uniform. For the latter discussions the arguments of the potential matrix V are omitted. Many relevant 
features of spatial evolutionary games can be investigated where the equivalent players are located on the sites 
of lattice with periodic boundary conditions. The investigation of these systems can be efficiently supported 
by the methods developed for studying models of solid state physics where all the relevant symmetries are 
satisfied. 

The latter simplifications exclude the analysis of systems with players interacting via non-symmetric 
games. At the same time, the systematic investigation of the multi-agent systems with non-symmetric 
games can be performed on bipartite networks that give us a convenient connectivity structure. For these 
multi-agent systems we distinguish two types of players ( e.g. males and females, buyers and sellers, etc.) 
who are distributed on a bipartite graph providing that each player plays games with players of opposite 
type. For these systems the connectivity network is divided into two subgraphs with sites x £ X and y £Y. 
The bipartite connectivity graph ensures that players at sites x € X interact exclusively with those at sites 
y £ Y, and vice versa. 

Such a situation can even be realized in lattice systems. For example, if the players are located at the 
sites of a square lattice, then the sites are divided into two equivalent sublattices corresponding to the white 
and dark squares of the chessboard. These systems will be discussed later in Sec. 16.II and IjUTl 

The first investigations of the multi-agent evolutionary games were performed for well-mixed populations 
corresponding to a complete connectivity grap h where interactions exist between all the possible pairs. Since 
the pioneering work of Nowak and May 


1211 ] the games are studied progressively on different lattices. 
For these spatial systems the sites are equivalent if periodic boundary conditions are applied for a finite 
square box of L x L = N sites and the interactions are restricted to the same size of neighborhood for each 
player ( e.g ., nearest or nearest and next-nearest neighbors). The main advantage of the lattice structured 
population is based on the translation symmetry providing good conditions for studying the consequences 
of some characteristic features in the structure of neighborhood. 

For real social systems, however, the irregular networks give a more adequate background for the investi¬ 
gation of multi-agent systems. This is the reason why the analysis of lattice systems was extended to different 
graphs representing diluted lattices 0SH real social networks S3] random graphs (8^, SH SH IH, [87] in¬ 


cluding different small-world 


and scale-free networks 
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4-1. Pure Nash equilibia in multi-player games 


In this section we briefly survey the relevant variations in the number of Nash equilibria when an TV- 
player game on a lattice (or network) is built from pair interactions between the neighbors. Regardless of the 
connectivity structure, the multi-agent systems have 2 N microscopic states (s) if all the pair interactions are 
characterized by a 2 x 2 game. Furthermore, if the possible transitions between two microscopic states are 
limited to those s — > s' where only one player modifies her strategy ( e.g ., s x —> s^) then the dynamical graph 
can be represented by an TV-dimensional hypercube 41, 92]. Such a dynamical graph is shown (for TV = 4) in 
Fig. [TT] where the nodes and edges of the four-dimensional hypercube are projected onto a two-dimensional 
plane from a suitable direction in such a way that the number of players playing the first strategy differs by 
1 from level to level. Notice that the parallel edges represent strategy changes for the same player while the 



Figure 11: The flow graph in a four player potential game if the potential is composed of two-player coordination games when 
players benefit from choosing the same strategy. 


strategy profiles are different on the neighborhood. If the strategy profile evolves along the directed edges 
then the system evolves into one of the two Nash equilibria. Such a behavior occurs when randomly chosen 
players are allowed to increase their income by choosing another strategy in the knowledge of their other 
possible payoff(s). The iteration of this process ends in one of the Nash equilibria that are absorbing points. 
Evidently the final state depends on the initial state. The latter method can be considered as a way to find 
pure Nash equilibria 93|, [44, 94, |29j. 

This method can be utilized for all the potential (and ordinal potential) games to find pure Nash equi¬ 
libria due to the absence of directed loops in their flow graph. The TV-player systems have only one Nash 
equilibrium if the interaction is composed of pair interactions belonging uniformly to either the harmony or 
prisoner’s dilemma games. 

In the example illustrated in Fig.[lT]there are only two Nash equilibria and each has a basin of attraction 
in the space of strategy profiles. A similar behavior occurs for larger number of players if all the possible 
pairs interact via a coordination game. 

More sophisticated processes occur in a spatial system where the players are located on a square lattice 
and play stag hunt games with their four nearest neighbors. In that case if we allow a randomly chosen 
player to increase her payoff by choosing the opposite strategy then the iteration of these steps drives the 
system into the homogeneous D state ( s x = D \/x) if the sufficiently large system is started from a random 
initial strategy profile for S < T — 1 and T < 1 (see the left plots in Fig. [12]). This is the preferred Nash 
equilibrium where U (s) reaches its maximum. This system, however, has many other Nash equilibria that 
may be achieved by the application of the above evolutionary dynamics if we choose different initial states as 
it is illustrated in Fig. [12] The reader can easily check that in the homogeneous C state (s,*, = C Mx) nobody 
can benefit from changing her strategy unilaterally. This feature implies that the latter Nash equilibrium 
has a finite basin of attraction if the formation of a D — D pair on two neighboring sites is not favored 
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Figure 12: Final spatial distributions of strategies (Nash equilibria) (lower snapshots) if the evolution is started from three 
different initial states (upper snapshots) when the players play stag hunt games (T = 0.5, S=-0.6) with the four neighbors on 
the square lattice. In the left case the system is started from a random initial state where D and C strategies are present with 
the same probability. In the initial state of the middle case Ds are created randomly with a probability of 0.02. For the right 
case the latter initial state is modified by adding two rectangular blocks of Ds. 


either. Consequently, the second Nash equilibrium can be achieved from such a random initial state where 
the ratio of D and C strategies is small and typically only isolated D strategies or small cluster of Ds are 
present in the sea of Cs. Such a situation is illustrated by the middle pair of snapshots of Fig. [T2] 

Despite the small ratio of D strategies, larger clusters of Ds may also appear in the random initial state, 
particularly if the system is sufficiently large [for a more quantitative analysis we suggest consulting surveys 
of percolation theory 95, [96|]. If the size of a rectangular cluster of Ds exceeds a threshold value, then this 
colony can remain alive or it can even expand if this process is supported by the presence of solitary Ds 
along the periphery as demonstrated in the plots at the right of Fig. [T2] Consequently, the selected Nash 
equilibrium depends on the initial state and also on the random strategy refreshment of players. 

Similar phenomena can be observed for the spatial stag hunt games if S' > T — 1. The only difference 
is that the roles of the D and C strategies are exchanged in the evolutionary processes. On the contrary, a 
basically different behavior occurs along the boundary S = T — 1 separating the above discussed regions of 
parameters where the homogeneous states form a frozen poly-domain structure as demonstrated in Fig. 1131 
As this is the typical behavior in the whole region of the hawk-dove games, we will not discuss the resulting 



Figure 13: Typical distribution of D and C strategies on a square lattice in a Nash equilibrium occurring in the stag hunt 
region for S = —0.5 and T = 0.5. 


processes separately. 

Within the region of hawk-dove games (T > 1 and S > 0) the lattice system has two equivalent sublattice 
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ordered states where the potential reaches its maximum. The analysis of these states requires the introduc¬ 
tion of two sublattices (X and Y) as described above. In the first ordered structure the sites of sublattice 
X are occupied by Ds , and C s are in the other sublattice. For the second ordered structure the sublattice 
occupancies are reversed. These long range ordered structures are equivalent (because of the translation 
symmetry) and can be observed in finite systems if the linear size L is even for periodic boundary conditions. 
If such a system is started from a random initial strategy distribution, then the random strategy update 
of players favoring the increase of their own payoff results in a frozen poly-domain structure as illustrated 
in Fig. [TT] Similar patterns are reported by Sysi-Aho et al. [97], who studied additionally the cases where 
second neighbor interactions are also taken into consideration. 




Figure 14: Typical strategy distributions on a 50 X 50 box of sites for Nash equilibria in a HD game on a square lattice for 
nearest neighbor interactions if S = 0.2 (left panel) and S = 0.4 (right panel) for T = 1.3. 


In the frozen poly-domain structures both types of ordered strategy arrangements are recognizable within 
the domains. Notice the absence of point defects inside the ordered spatial regions where the players are 
satisfied and will not modify their strategy. The opposite ordered structures are separated by boundary 
layers. The average composition of these boundary layers, however, depends on the sign of S — T + 1 as 
demonstrated by two opposite examples in Fig. [IT] 

Evidently, the above-mentioned difference in the boundary layers vanishes if 5 — T + 1 = 0 because the 
interface becomes a mixture of two types illustrated in Fig. [II] In this special case there are many players 
along the interfaces who have the same payoff for the opposite strategies. If we additionally allow these 
players also to modify their strategy, then these additional stochastic events in the strategy update yield 
a domain growing process that drives the system into one of the homogeneous sublattice ordered strategy 
arrangements. This evolutionary process is demonstrated in Fig. [15] where the time is measured in the unit 
of Monte Carlo steps (MCS) [within one MCS each player receives a chance once on the average to change 
her strategy]. 



Figure 15: Three snapshots at times t = 10, 40, and 160 MCS (from left to right) illustrating the evolution of strategy 
distribution on a 80 X 80 box of a larger system for S = 0.4 and T = 1.4. 


Finally we mention that the latter additional options result in a similar domain growing process in the 
stag hunt region, too, if S — T + 1 = 0. These ordering processes will be discussed later when studying 
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the consequences of a more general class of stochastic evolutionary rules. Finally we mention that similar 
phenomena occur on bipartite networks. 


5. EVOLUTIONARY POTENTIAL GAMES 


Evolutionary games are devoted to explore the consequences of different evolutionary rules defining the 
ways how individuals can modify their own strategies in order to have a higher payoff. The first efforts were 
mainly focused on the maintenance of cooperative behavior in the social dilemmas where the results strongly 
depend on the set of strategies, the connectivity structure of players, and the dynamical rules that may 
involve variations in all ingredients of the mathematical models (including the current strategy, interaction 
or learning networks, dynamical rules, personal features, and their combinations). As mentioned in the 
Introduction the concept of evolutionary games was originally suggested and applied by Maynard Smith 
and Price 11], Maynard Smith Q, Axelrod and Hamilton fuj . Axelrod [16], Hofbauer and Sigmund Q to 
provide a mathematical background for the description of the Darwinian evolution for biological and social 
systems. This is the reason why most of the investigations during the last decades used dynamical rules 
based on the adoption/imitation of the more successful strategy (or biological species). As a result different 
dynamical rules are introduced and studied systematically for a wide scale of conditions [for a survey see 

BElil]. 

For example, in lattice systems the ’’imitating the best” rule [20] dictates the deterministic adoption 
of that player’s strategy who received the higher payoff in the neighborhood; the ’’imitating the better” 
rules allow to adopt any better strategy from the local neighborhood with a probability dependent on the 
iH, 18] ■ B eside s it there are other imitation rules [introduced and studied e.g., by Nowak 
; Szabo and Toke fiool : Alonso-Sanz et al. (101| : Masuda and Aihara j84j; Ohtsuki and Nowak 


neighbors’s payoff 
et al. 


[102] , Wild et al. [103 ]; Wu et al. [104 ]; and many others] when the players can adopt the strategy even from 
a worse neighbor though the more successful strategies continue to be preferred. 

For all types of imitation rules the system has absorbing states because this rule cannot recreate a 
strategy that has become extinct previously. For example, if one of the homogeneous stat es is a ppro ached 
then the system remains there forever. On the contrary, the introduction of mutation [see jl05| . [106], and 


107]] is capable of reviving strategies that have vanished previously. 


For the rules based on fictitious game (as detailed in the SectionQj]) each player is capable of determining 
her own payoff for her strategies with an assumption that their co-player will not change their strategies. For 
the ” best response” version of these rules the player chooses those strategies that provide the highest payoff 
under the mentioned conditions. If more than one ’’best response” exists, then the player selects one of them 
at random. The smoothed versions of this rule allow the player to select all her possible strategies with a 
probability dependent on the payoff variation. The concept of potential games becomes fruitful for those 
special versions of these dynamical rules which drive the system into the so-called Gibbs st ate ch aracterized 
by the Boltzmann distribution [for an English translation of Boltzmann’s original work see Il08] ]. For these 
special cases we can utilize the enormous amount of results obtained in statistical and solid state physics to 
explain phenomena emerging in social and biological systems, too. 

Evidently, the latter evolutionary processes based on fictitious games are relevant in human systems. 
Similar situations can occur in biological systems if the species creates a large number of mutants with a 
survival probability dependent on the neighborhood. We have to emphasize, however, that the quantitative 
analysis of the social and biological systems requires relevant extensions of the possible rules opening new 
challenges for both the non-equilibrium statistical physics and super-statistics. 


5.1. Evolutionary rules leading to Boltzmann distribution 

In the field of game theory Blume [28] described first that for the application of the so-called log- 
linear or logit evolutionary rules the potential games evolve into the Gibbs state assuring the validity of 
the direct application of the concepts and tools of equilibrium statistical physics. The basic idea, as an 
efficient Monte Carlo algorithm to determine average values for the Boltzmann distribution in many particle 
system, however, was developed by a group conducted by Metropolis at Los Alamos in the 1950s. They 
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recognized that instead of evaluating the Boltzmann factor for a randomly chosen microscopic state it is 
more efficient to use a Markov chain algorithm in which a state will be chosen with a probability defined by 
the Boltzmann distribution. No wada ys this method is named importance sampling. A similar evolutionary 
rule was suggested by Glauber [ 109| who wished to study time-dependent processes for the Ising model 
whose dynamics is not prescribed by the quantum mechanics. Since that time several other algorithms 
were suggested and used widely for the investigation of equilibrium and non-equilibrium phenomena in 
many-particle systems. 

Now we describe the above-mentioned rules in historical order by using the terminology of evolutionary 
game theory. The system is started from an arbitrary (typically random) initial state and during the 
evolution the following elementary steps are repeated. A randomly chosen player x can modify her strategy 
s x to s(, chosen at random from her options with a probability w(s x — > s(.) dependent on the payoff variation 
A u = u x (s^.,s_ x ) - u x (s x , s_ x ) = V(s' xl s_ x ) - V (s x , s_ x ). 

For the Metropolis algorithm the player always adopts the new strategy s' x if it is worthy (A u > 0). For 
the opposite cases the player accepts the strategy with a probability decreasing exponentially with the 
possible loss, that is, 


w 


(M) 


(s x s x ) 


1 for A u > 0 , 
e A u/k fo r Am < 0 


(69) 


where K denotes the magnitude of noise (or temperature in physical systems). Notice that this rule agrees 
with those we used in Sec. 14.11 if K —»• 0. 

For the Glauber dynamics 


w {G) (s s' x ) 


1 

l + e~ Au /K' 


(70) 


This evolutionary rule is convenient for the analytical calculation as demonstrated in 


The transition rate for the logit rule 2 


ne analy 
AS. 110 


is defined as 




( s * "4 s'J — A a 


d s A s -x)/k 


E s 


j(s",s - X )/K ' 


(71) 


where the summation runs over all the possible strategies of the player x and 0 < A x < 1. The choice of 
A x = 1 provides the most efficient algorithm in the Monte Carlo simulations. 

The most relevant common feature of the above dynamical rules is that they all drive the potential games 
into the Gibbs state where the probability p( s) of a strategy profile s is given by the Boltzmann distribution 
as 

P(s) = (72) 

where 

Z = ^2e u(s)/K . (73) 

S 

The normalization factor Z is called the partition function and plays a crucial role in statistical physics 
because the partial derivatives of In (Z) reproduce relevant thermodynamical quantities. For example, the 
average value of the potential can be given as 


9 In (Z) 

~ d(l/I<) ■ 


(74) 


lllj . Baxter 


The relevance of these features is enhanced by the exact solutions obtained by Eggarter 
E2, Yang |l 13j] for several connectivity structures. 

In addition, the conditions of detailed balance are also satisfied between all the possible forward-backward 
transition pairs. As a result, for the Boltzmann distribution the transition s x — > s' x and its backward version 
(s(, —> s x ) appear with the same frequency, i.e., 











for V x, V s x , s' x , and V s_ x . These conditions are satisfied if 


w(s x —> S^S-a;) _ [u(b')-U(b)]/K 

w(s' x -> s x ,s_ x ) 


( 76 ) 


Due to this relationship the final stationary state (the Boltzmann distribution) and also the satisfaction of 
detailed balance remain unchanged if the transition rates w(s x —> sj.) and w(s' x —► s^) are multiplied by the 
same coefficient that can be varied from pair to pair of strategy profiles involved [see X x in Eq. CD]. 

For all the above dynamical rules only one player has modified her strategy in the consecutive elementary 
steps. For the Kawasaki dynamics [114] a randomly chosen player x exchanges strategy with one of her 
neighbors with a probability dependent on their summarized payoff variation that is not affected by the 
interaction between them. The payoff variation arising from the rest of interactions can be built up from 
two consecutive individual strategy changes as before. Consequently, the summarized payoff variation is 
identical to the potential variation A U = U(s x , s y , S- X , y ) — U(s y , s x , S— X<y ) and the transition rate can be 
given by the following expression: 


w 


(K), 


s s *) — 


1 


1 + e ~ AU / K 


(77) 


used in many solid state applications when studying diffusion or transport in systems where the number of 
particles (strategies) is conserved. For the two-strategy systems the corresponding dynamical graph is given 
by the nodes of an TV-dimensional hypercube and the edges are defined by those diagonals of faces where 
two neighboring players exchange their different strategies with each other. 


5.2. Statistical physics and thermodynamics 

For the microscopic description of a multi-agent system we should define the individual strategy for all 
participants, which requires a huge amount of data. At the same time, for the macroscopic description 
of the same system we use only a few parameters, e.g. average portion of a strategy, average payoff and 
noise level, as it happens in the thermodynamic description of a gas when using the concept of pressure, 
temperature, and density, without any knowledge about the position and velocity of atoms contained in the 
gas. Statistical physics gives a general framework to describe the relationships between the macroscopic 
(thermodynamic) quantities that are influenced by the microscopic interactions. 

The fundamental assumption of the equal a priori probability of the accessible microscopic states serves 
a basis for the statistical description. According to this approach, those macroscopic states can be observed 
in the stationary state that are realized by the largest number of microscopic states. Due to the law of 
large numbers the latter approach defines well the macroscopic state in the limit N —> oo as d etailed in the 
textbooks of statistical physics [see e.g. EE EH] based on early works by Boltzmann E3; Gibbs }l!8| 
and Szilard EE- 

Now we briefly survey the relevant mathematical background followin g a moder n formalism suggested 
by Jaynes [l20l ] and using the terminology of evolutionary game theory 28L 121 . 122 ] as before. Within this 
framework the system is described by the p(s) probability distribution of the micr oscop ic states s and the 
problem of prediction becomes equivalent to the maximization of Shannon entropy 123] 


S = -^2 p(s) In p{s) 


(78) 


under several constraints. The crucial problem of extending the maximum entropy principle to non-physical 
systems lies in the adequate choice of constraints. This approach is used successfully in the analysis of 
complex 124| and chaotic systems. 

In the present case p{ s) is normalized, that is, 


J2p( s ) = 1 


(79) 
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( 80 ) 


and we can fix the average value of the potential as 

J2p( s )U(s) = U. 


Using the standard method the extremum under the above constraints can be evaluated by solving the 
following set of equations: 


d 


dp(s) 


S + Ai p ( s ) + A 2 p( s ) u ( s ) 


= 0 


(81) 


where the Lagrange multipliers Ai and A 2 are determined by taking into consideration the conditions G9 
and m- Straightforward calculations yield that the optimal probability distribution becomes equivalent to 
the Boltzmann distribution for A 2 = 1/K [see (l72l) and ([73)1 ] and for a suitable choice of the average value 
of the potential (U) depending on K (strictly) monotonously. 

In the real utiliz atio n of the above extremum principle, it is more convenient to perform a Legendre 
transform 125 . Hii lml. Using the analogy of the Helmholtz free energy in thermodynamics, we introduce 
now a thermodynamic potential 

<f> = U + KS (82) 

that has a maximum in the equilibrium state for a fixed I\ characterizing the temperature in statistical 
mechanics or noise level in evolutionary potential games in the sense defined by the logit dynamical rules 
[see e.g., (I7T1) ]. 

Notice that for spatial systems with short range interactions both the entropy S and the average value 
of potential (as well as <f>) are proportional to the system size N and are considered as extensive quantities 
in the thermodynamical descriptions. 

Counterexamples are systems with long-range interactions or systems on some types of small-world 
networks where we cannot apply the traditional methods directly. In the wide field of statistical physics 
there are some new directio ns an d approaches to describe the macrosco pic behavior of the mentioned systems 
when considering complex 124 1. or other dynamical systems 128l 123 , and those non-equilibrium systems 
that behave like the superpositions of the Boltzmann distributions 130), 131]. Surveying the latter directions 
goes beyond the scope of the present work, our attention will be focused on the description of phenomena 
occurring in large spatial systems. 

In the following sections we will draw a parallel between the kinetic Ising/Potts models and social systems 
for a particular evolutionary rule. In fact, for both types of these systems the applicability of thermodynamics 
is limited to the thermostatics. According to quantum mechanics each microscopic state of the Ising type 
systems is an eigenstate and remains unchanged forever. The time-dependence in the kinetic Ising model 
is introduced by assuming a uniform interaction between the spins and a heat reservoir characterized by 
the temperature K. Under these conditions the total energy of the Ising system fluctuates and the first 
and second laws of thermodynamics become meaningless. The failure of the third law of thermodynamics 
is related to the high degeneracy of the ground state as discussed later for several cases. 

In social and biological systems the application of a similar stochastic evolutionary rule is justified by 
its mathematical simplicity and the analogies we can exploit in the interpretation of the phenomena. In 
these systems the value of K can characterize: (i) the fluctuation of payoffs; (ii) the errors made by players 
during their decision process; and (iii) the magnitude of risk the player accepts in the hope of finding a 
better (long-term) solution. In social systems the value of K can be even considered as a personal feature 
and also the subject of the coevolutionary process. The first simulations in such coevolutionary models have 
indicated the homogenization of K for the coexistence of strategies within the prisoners dilemma region if 
both the strategy and K adoption are controlled by pairwise imitation rule |l32| . Similar results for Glauber 
type rules would increase the relevance of thermostatics for these systems. 

5.3. Consequences of the extremum principles 

The principle of the maximum entropy serves as a mathematical background describing the intimate 
relationship between statistical mechanics/physics and thermodynamics. This principle explains the laws 
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of thermodynamics as well as the relevance of the Boltzmann distribution and connects the meaning of the 
Lagrange multipliers and the intensive quantities of thermodynamics. The Boltzmann distribution itself 
implies several general relationships among the first and second partial derivatives of the thermodynamic 
potentials 125 . 1 33l . 1 .'1 11 or the partition function as illustrated by Ecp (l74l) . Examples are the Gibbs-Duhem 
relations and the Gibbs’ phase rule. This knowledge is summarized in equation of states quantifying the 
relations among the intensive and extensive thermodynamical quantities for macroscopic systems composed 
of different atoms, ions, and molecules. 

The simplicity of the linear response th eory and the fluctuation-dissipation theorem are also direct 
consequences of the Boltzmann distribution 135| where the effect of a small perturbation obeys a simple 
expression in linear approximation. 

In addition to the general thermodynamical relationships, the extremum principles can also be used 
to evaluate the thermodynamical quantities in the knowledge of the microscopic interactions. We briefly 
survey now the essence of the cluster variation methods providing a general framework for the traditional 
mean-field and pair approximations in the approximative description of the lattice systems. For the sake of 
simplicity our description is focused on systems with interactions between the equivalent neighboring players 
located on the sites of a lattice. Within this approach the translation invariant microscopic states can be 
described by configuration probabilities on a small cluster of neighboring sites. For example, pi(si) defines 
the probability of the strategy si at each site of the lattice while P 2 (si, S 2 ) describes the probability of finding 
si and S 2 strategies on two neighboring sites. For rotational invariant arrangement of strategies ( e.g ., on 
square lattices where the horizontal or vertical directions are equivalent) we can assume that p 2 (si,S 2 ) is 
independent of both the position and direction of the pair. These quantities are normalized, that is, 


E*(-o = 

Si 

E P2(si,S 2 ) = 1 , 


(83) 


and satisfy the compatibility conditions: 

J2p2( Si ,S2) = EP 2 ( S2 ’ S l) = Pl( s l) ■ 


(84) 


Evidently, one can use larger clusters of n sites to describe the stationary state in this system and the 
corresponding quantities satisfy similar compatibility conditions. The larger the cluster we study, the more 
accurate is the approximate solution [1361 11371 ]. 

With the application of Bayes’ theorem we can build up the configuration probabilities for a large cluster 
as a product of configuration probabilities of smaller clusters as detailed in the paper by Gutowitz et al. 

138j . This approach allows us to derive adequate approximations for some other quantities, e.g., correlation 
function and correlation length 22j. Furthermore, the application of Bayes’ theorem plays a fundamental 
role within the dynamical cluster techniques (developed to study stationary states in non-equilibrium lattice 
systems) when a set of coupled differential equations is derived by taking into consideration the contribution 
of all the elementary steps |l39( | . In comparison to the dynamical cluster techniques the so-called cluster 
variation method provides a more convenient way to evaluate the configuration probabilities. With the use 
of these quantities one can express both the entropy and the average value of the potential for a given lattice 
structure. For example, on a d-dimensional cubic lattice the value of U is expressed as 


U = Nd E Vs 1 s 2 P2(s 1 ,S2) , 

<Sl,S2 

while the entropy can be approximated at the level of pair approximation [M Il4ll | as 


S ~ <S (2p) = -Nd E P2(si,s 2 )lnp 2 (si,S2) + N(2d- 1) Epi( s i) ln Pi( 5 i) • 

si,s2 si 


(85) 


( 86 ) 
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When using larger cluster sizes one can give better approximations for the entropy 142l . 143], while ex¬ 
pression (1851) remains unchanged for the case of pair interactions. For the cluster variation methods the 
thermodynamic potential [ e.g ., $ as defined by m] is expressed as a function of configuration probabili¬ 
ties and its maximum value is determined by varying these parameters. The standard calculation may be 
simplified by considering only the independent parameters of the configuration probabilities. For example, 
if only two strategies are allowed, si, S 2 = 1,2, then the one- and two-site configuration probabilities can be 
expressed with only two parameters as 


Pi(l) = c, 

Pi (2) = 1 - c, 

P 2 (i,l) = q, 

p 2 (l,2) = c-q, (87) 

4 * 2 ( 2 , 1 ) = c — q, 

P 2 ( 2,2) = 1-2 c+q 


where c can be interpreted as the portion (or density) of strategy 1 in the whole system and c — q describes 
the density of domain walls separating homogeneous territories. In that case the equilibrium values of c and 
q are given by the solution of the following equations: 


and 


d$(c,q) 

dc 

d${c,g) 

dq 


= 0, 


= 0, 


( 88 ) 


(89) 


where <I>(c, q) = U(c,q) + KS(c,q) is obtained by substituting Eqs. (IH71) into (1M1) and (l85]l . In general, 
Eqs. (HU) and (l89l) have more than one possible solutions (0 < 4 * 2 ( 52 , 52 ) < 1)- In the latter case the real 
solution is the one where <h(c, q) reaches its maximum. 

Thus the cluster variation method gives us a unique way to evaluate the equilibrium value of c and q 
together with all the related quantities as a function of K. Besides it, one can derive approximate phase 
diagrams if the actual model has several solutions with differe nt sy mmetries, as it happens frequently in solid 
state systems [examples and further references are given in 
the average total payoff can be given as 


144 Il45j ]. In evolutionary games on lattices 


A = Nd^ A slS2 p 2 (si, s 2 ), 

Sl,S 2 


(90) 


in the knowledge of the payoff matrix A and the nearest neighbor strategy configuration probabilities. Apart 
from the symmetric case A = A T , the quantity A has no analogous concept in solid state physics because 
it may contain terms neglected in the evaluation of the pair potential V. At the same time, this quantity 
plays a key role in the investigation of social dilemmas. 

It is worth mentioning that the cluster variation method at the level of pair approximation reproduces 
the exact result for the one-dimensional systems with nearest neighbor interactions. For these systems the 
probability distribution of strategies can be given by the Bayesian formula as 

x=N—l / \ 

p(s) =p 2 (si,s 2 ) n P2 Sx l Sx , +1 , (91) 

Pl{s x ) 

s x ,x =2 v ' 


that satisfies the compatibility conditions, that is, any pair configuration probability can be reproduced by 
summing over the rest of sites, and the resultant entropy becomes equivalent to (1861) for d = 1 in the limit 
N —> 00 . 
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Equations from (1551) to (15T1) give a mathematical basis for the use of the cluster variation method at the 
level of two-site approximations. This calculation becomes simpler at the level of one-site approximation 
(that is equivalent to the traditional mean-field approximation), when we have only one parameter (c) to 
be determined according to (188|) . The corresponding expressions for U and S can be derived with the 
assumption p 2 (si,S 2 ) = Pi(si)pi(s 2 ) that simplifies the entropy as 

S ~ <S (lp) = -AT^pi(si)lnpi(si). (92) 

Sl 

It is emphasized that the one-site approximation ignores the role of the topological features of the connec¬ 
tivity structure. More precisely, only the number of co-players is involved in U. This approach may give an 
adequate description of phenomena for those (homogeneous) connectivity structures where the number of 
neighbors is large enough. For the square lattice the results of the one- and two-site approximations will be 
contrasted with Onsager’s exact result in the following section. 

The generalization of the cluster variation method for larger number of strategies and/or for larger size of 
clusters is straightforward. When increasing these parameters one can study a rich variety of sophisticated 
phenomena although the calculations become more complex and time-consuming. At the same time we 
can reduce the number of independent parameters by taking into consideration the additional symmetries. 
From this point of view the use of the corresponding three- and four-site approximations is advantageous 
on triangular and square lattices. 

Finally we mention that the pair approximation can be applied successfully on Bethe lattices where this 
approach gives a more accurate prediction as illustrated by Vukov et al. [86| who compared the analytical 
approximate results with Monte Carlo simulations performed on a locally similar random regular graph for 
large sizes. 

In the above description we have assumed the equivalence between both the players and their location. 
There exist, however, several spatial structures ( e.g. non-Bravais lattices) or games (e.g., matching pennies) 
where we should distinguish the sites and/or players. The cluster variation method may be extended to 
these cases by introducing the configuration probabilities on several types of one-, two-, or n-site clusters 
by considering also the sublattice structure. 


6. ISING MODELS 


The inv estig ation and application of the Ising type models have a long history as su rvey ed by Brush 
tl4fij | . Niss [l47l ll48L ll49j |. Sornette (TEcj. The original idea was first described by Lenz 151 1 as a simple 
model to study the ferromagnetism. The one-dimensional version of this lattice model was studied by Ising 
Il52| . who was the PhD student of Lenz. Unfortunately, the one-dimensional model is not suitable to describe 
the ferromagnetic-paramagnetic transition observed in magnetic material when increasing the temperature 
EH EH because the domain walls, that c an b e considered as point defects in the one-dimensional systems, 


prevent the formation of long range order [1151 . 


The exis tence of the phase transition was shown by Peierls 155| using an a rgument later improved 


by Griffiths [l56| . The approximate methods (like mean-field approximation 157 and pair approximation 
140|) have confirmed the presence of a continuous phase transition. The exact solution in the absence of 
exter nal magnetic field on the square lattice with nearest neighbor interactions was obtained by Onsager 
[l58j . In the following decades this magnetic model was extensively investigated on different lattice structures 
assuming ferromagnetic or anti-ferromagnetic interactions. Due to its simplicity and the knowledge of the 
exac t solution the model was frequently used to check the accuracy of different approximation techniques 
EH- Since that time the Ising models have been fundamental mathematical mo dels in statistical physics 
representing a robust universality class of critical phase transitions EH EH, Ib il Ttiii 


The original Ising model can be directly applied to explain the magnetic behavior only for a few materials. 
On the other hand, the equivalent lattice gas models are widely used to derive theoretical phase diagrams for 
alloys [l63j . met al-hydrogen systems 164 1, solid electrolytes |165| . intercalation fl66l | or non-stoichiometric 
compounds [167{ , etc. The analogy to the Ising models is due to the fact that in all these models only two 
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states of the lattice sites are distinguished. For example, in the Cu-Au alloys either Cu or Au atoms can be 
at one site x; in the palladium-hydrogen systems an interstitial void x of the metal matrix can be empty or 
occupied by a H atom; in the superionic conductor Agl, iodide ions form a rigid body-centered cubic lattice 
and the smaller mobile silver ions can be present or missing in a tetrahedral void represented by the site x. 
The total energy of these systems can be built up from pair interactions between the neighboring sites in 
the knowledge of the physical properties. 

The flexible interpretation of the Ising model implies its applicability to many other systems involving 
high-energy physics [for a survey see t he review by Pelissetto and Vicari jlfifij]] and social models with 
players located on a network jlfifA I170j . The latter situation occurs when the connected individuals can 
choose between two options, e.g., using meter or yard as unit of length; following drive left or drive right 
rule in traffic; using Windows or Linux operating systems, etc. All t hese example s represen t the so-called 
coordination game and are analogous to ferromagnetic systems jl69l . ll7lL Il72l Il73l . Il74 . 175l 


176]. 


6.1. Systems equivalent to Ising models 

In the magnetic Ising model the spin variables a x = +1,-1 refer to upward and downward magnetic 
moments at site a" of a network. The strength of the spin-spin interaction between the neighboring sites x 
and y is denoted by J xy and we can assume the presence of a site-dependent external magnetic field h x . For 
any microscopic state a = (<r i,... , er/v) the Hamiltonian (or potential energy) function of this multi-spin 
system is given as 

H (< 7 ) — ^ ] JxyGx&y ^ ] h x (X x ( 93 ) 

(x,y) x 

where the summation runs over all nearest neighbor pairs (denoted by (x , y) as it is done for the multi-player 
games). For ferromagnetic interactions the so-called coupling constants are positive ( J xy > 0) and the system 
achieves one of the ordered ground states with minimal energy in the ferromagnetic phase (a x = +1 at each 
site if h x > 0 or a x = — 1 if h x < 0. A similar optimal strategy distribution occurs in a multi-agent system 
with coordination type pair interactions if one of the two strategies is preferred by suitable self-dependent 
payoffs. 

In order to clarify the relationship between the Ising models and the multi-agent, two-strategy, potential 
games first the site energies are shared among the pair interactions as 

h x = Y J KW), (94) 

y' 


where the summation runs over y' which are the interacting neighbors of x. In that case the Hamiltonian is 
composed of pair interactions involving the corresponding contributions of the site energies from both sites 
in the following way 

H(<x) = H xy( ax ’ a y) (95) 

(x,y) 


where 


Hxy {eT :r . (7y ) Jxy(X x (Xy U ;; ■ 


that may also be written as 


h'y{x)<Jy 


(96) 


Hxy {(X x ; &y ) — Jx 


,G X (Ty 


K(y) + K ( x ) 


{(X x T ®y ) 


K(v) - K( x ) 


( a X (Xy). 


(97) 


The latter form of the pair interactions is convenient for the comparison with the potential of a nonsymmetric 
2x2 game given by Eq. (1M1) because the first term is analogous to the coordination type interaction as the 
four possible values of the product a x a y define a matrix equivalent to f(4). Similarly, in matrix notation 
the second (third) term of (11)71) becomes equivalent to the second (third) term of the expression (1641) with 
a suitable choice of the coefficients. Using this analogy a multi-agent, two-strategy, potential game can 
be mapped onto an Ising model on the same network with an adequate choice of the parameters J xy and 
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h' x (x') whereas the local site energy is defined by (TO1) . More than one game can be mapped onto the same 
generalized Ising model by varying the game parameters under the conditions of (1941) . 

In a human system it is natural to assume distinct payoff's for each pair, therefore in the corresponding 
Ising models both the coupling constants J xy and the local magnetic fields h x can be considered as random 
parameters. The resulting random Ising models or spin glasses will be discussed briefly in Sec. 17.101 

The traditional Ising model was introduced to study solid state phenomena where the interacting particles 
are identical and the applied magnetic field is homogeneous, consequently J xy = J and h x = h , and the 
connectivity network is a translation invariant lattice where each site has z neighbors. Similar conditions can 
be satisfied for multi-agent evolutionary games with equivalent interactions between the players residing on 
the sites of the same lattice. In that case the local site energy is shared equally among the z pair interactions 
becoming 

H xy i (Ty ) — J &x&y (*4c T ) ■ (98) 

Z 

For the quantitative relationship now we compare —H xy with the pair potential (1621) obtained in the notation 
(l59l) of social dilemmas. Drawing a parallel between these two systems, the a x = +1 ( a x = — 1) spin state 
corresponds to the strategy s x = D (s x = C ). For this convention one can deduce the following relationship 
between the parameters of the Ising model and social systems: 


J = 


1-T-S 

4 


and 


h 

z 


T-l-S 

4 


(99) 


that is, the Ising type spin-spin interaction is analogous to the coordination game whereas the magnetic 
field plays the role of the driving force validating the risk dominance. 

These linear relations determine the values of J and h as a function of the payoff parameters T and S. 
Figure [16] shows the orthogonal Descartes coordinate axes for the h/z and J (denoted by red dashed and 
blue dotted lines) on the S — T plane we used in Sect. 13.31 for the classification of games according to the 
flow graphs and Nash equilibria. Now the S — T (or h — J) plane is divided into three regions with respect 



Figure 16: (Color online) Mapping the Ising parameters J and h onto the S — T plane of payoffs in the notation of social 
dilemmas defined by Eq. J59J. 

to the thermodynamically stable states of the Ising model in the low noise limit. Within the white territory 
a ferromagnetic, i.e., homogeneous spin down (a x = —1), is stabilized that refers to uniform cooperation 
in games (where T < 1 and S > T — 1) involving the harmony games and the upper half of the stag hunt 
games. The opposite ferromagnetic state ( a x = +1 or s x = D) occurs within the dark region (S < 0 and 
T > S + 1) representing the second half of the stag-hunt games and the prisoner’s dilemmas. 

Within the third region (T > 1 and S > 0) an anti-ferromagnetic spin arrangement is realized when 
the neighboring spins point to opposite directions. On the square lattice the up and down spins form a 
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chessboard like pattern which we used in Fig. [TH] for the illustration of this spatial structure. This anti¬ 
ferromagnetic spin arrangement is twofold degenerated as the simultaneous reversal of all spins does not 
affect the total energy (lil3l) even in the presence of h. On the other hand, the application of a magnetic field 
h influences the local stability of spin states and this is the reason why a ferromagnetic phase is stabilized 
if \h\ exceeds a threshold value when J < 0, as it is also discussed for games |79j. These thermodynamically 
stable (ordered) states are realized by the logit rules in the limit K —» 0 and coincide with the preferred 
pure Nash equilibria determined by the strategy pair (i,j) for which Vij is maximal. 

In the present system the highest average (total) income can be achieved by uniform cooperation if 
2 R > T + S, otherwise the sublattice ordered arrangement of the C and D strategies provides the highest 
total payoff. Despite this total payoff optimum the system falls into the state of ” tragedy of the community” 
within the dark territory indicated in Fig. 1161 In the zero noise limit the latter structure occurs for h < 0, 
except the sublattice ordered region where the sufficiently strong repulsive interactions prevent the formation 
of a homogeneous state. In some sense the negative magnetic field can be interpreted as a force driving the 
social system into the state of the tragedy of the community. 

Notice furthermore that the case of J = (1 — S— T)/4 = 0 refers to the absence of interactions between the 
neighboring spins. In such situations the independent spin reversals are controlled only by the magnetic field 
h and the system is equivalent to the collection of non-interacting spins that can be described analytically 


when considering onlva single spin. The corresponding games are called ”equal-gains-from-switching” [177 


45| in the literature. The donation game [9| represents one of the well investigated 


or ’’dummy” games 
examples. 

The analytical study of the anti-ferromagnetic spin arrangements requires the division of the lattice into 
two sublattices (x £ X and y £ Y) in a way discussed in Sec. SI Additionally, the concept of the Ising 
model was extended by introducing staggered magnetic fields that are uniform within a sublattice. Thus, 
h x = h + h s and h y = h — h s where h s defines the strength of the staggered magnetic fields. Now the pair 
interaction 

h, . h„ 


Hxy{&xi@y) — J (J, 


x u y 


- [px 

z 


+ a y) _ ipx a y )■ 


( 100 ) 


has a third term quantifying the effect of the staggered magnetic field that favors one of the sublattice 
ordered spin arrangements. The Ising model with staggered magnetic field is equivalent to an evolutionary 
potential game on the same bipartite network. Comparison of the matrix notation of Eq. (11001) with the 
2 x2 potential matrix (1M1) explains that the staggered magnetic field 


= a'(7), 


( 101 ) 


where z denotes the number of neighbors in the translation invariant connectivity structure. 

The introduction and application of the staggered magnetic field have no practical importance in the 
investigation of the magnetic materials. This quantity becomes relevant from theoretical point of view when 
justifying the equivalence of order-disorder transitions observed for the ferromagnetic and anti-ferromagnetic 
Ising models when varying the temperature. To be more precise, the effect of h s on the ferromagnetic state 
(J > 0 and h = 0) is equivalent to the effect of h on the anti-ferromagnetic state. 

Finally we mention that the sublattice-dependent local site energy plays an important role in those lattice 
gas models where two types of sites are distinguished. For example, in the crystals of Pd-H system [1641 ]. 
CaF 2 [l65| . and alkali-fullerides 145 1 the mobile atoms/ions can stay both in the tetrahedral and octahedral 
sites of the face-centered cubic lattice formed by the rigid components. 


6.2. Potts models 

The two-state Ising model without magnetic field was extended by Potts 178] who introduced a lattice 
model with n (n > 3) equivalent states at each site. Previously Ashkin and Teller [mj studied a particular 
four-component version of the Ising model on a square lattice. The lattice model with n s tates for a common 
interaction energy between the different local states was also suggested by Kihara et al. Il80( . The name of 
Potts model was proposed by Domb [l8l | and used worldwide since that time. For a comprehensive survey 
of the early results we can suggest consulting the review by Wu jl82 |. 
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In the terminology of game theory the Potts models correspond to multi-agent games with symmetric 
n-strategy pair interactions where the payoff matrix and potential are given bynxn unit matrices, that is, 
both the payoffs and potential matrices V can be expressed by Kronecker’s delta as 

Aj/j = Bij = Vij (x. y ) = JSij , (102) 

where i,j = 1,... ,n. For J > 0 this pair interaction is equivalent to a generalized coordination game with 
n equivalent Nash equilibria when the players choose the same option. This game involves all the d (p,q) 
[1 < p < q < n; see Eq. (1411) ] coordination type interactions with the same strength. 

In the literature of physics the Potts model is extended by introducing an external magnetic field (h > 0) 
or a self-dependent component preferring one of the strategies. The presence of this term suppresses the 
critical transition in a way as it occurs in the Ising model. 

The main features of the corresponding ferromagnetic Potts model are well described in the field of 
statistical physics. This set of models is used for the classification of the universal behaviors appearing in the 
order-disorder transitions when the temperature (noise) is increased. Due to its simplicity the investigation 
of the Potts model played an important role in the exploration of the critical phase transitions and it became 
a key model when testing different methods and approaches. It turned out that the critical behavior of the 
Potts model is richer and more general than that of the Ising model. 

Originally the Potts model was considered as the simplest mathematical model exhibiting an order- 
disorder phase transition for n > 2 while its kinetic versions allow the investigations of the time-dependent 
ordering processes. In the light of the robustness and universality of critical phase transitions it was later 
realized that in many substances the phase transitions can be interpreted by the application of various 
Potts models. For example, three-fold de generated ordered structures can be formed by atoms adsorbed on 
single crystal surfaces (for 1/3 coverage) 11831 Il84j or by mobile ions in layered solid electrolytes (e.g. silver 
/3-alumina) jl85| . A large variety of the experimental realizations of the two-dimensional, n = 4 Potts model 
was discussed by Dornany et al. [l86l |. 

In Sect. 17. 81 we will discuss briefly a lattice gas model for illustrating the emergence of equivalent ordered 
structures on a square lattice. The general and universal properties of the order-disorder phase transitions 
for the n-state Potts model will be detailed in Sec. 17.41 while the most relevant features of the ordering 
process are surveyed briefly in Sec. [8] 

7. FEATURES OF ISING MODELS 

Despite its simplicity the Ising model can be used to demonstrate a surprisingly wide range of phenomena 
including different types of ordered structures in the stationary state when varying the range and strength 
of interactions on lattices or graphs. Besides it, the kinetic Ising model is capable of illustrating a large 
number of dynamical processes how they evolve towards the final stationary state. Now we briefly survey 
the most relevant phenomena that can help us understand the behavior of the more general evolutionary 
games on graphs. 

7.1. Spontaneous symmetry breaking 

In the absence of a magnetic field h the ferromagnetic Ising model with nearest neighbor interactions 
on the square or cubic lattices has two equivalent ordered states that appear with the same probability 
according to the Boltzmann distribution d. Similarly, the Boltzmann distribution predicts the presence 
of the two equivalent anti-ferromagnetic ordered structures with the same probability even in the presence 
of a magnetic field if its strength does not exceed a threshold value. On the contrary, in the practice we see 
only one of the ordered structures in sufficiently large systems. 

To resolve the above discrepancy Fig. [TT] shows a typical time-dependence of the frequency of C strategy 
as a function of time in one of the two sublattices of a small square lattice for the hawk-dove region at a 
low noise level when the sublattice ordered strategy arrangements are favored. Figure [l7] illustrates that 
the system alternates between two ’’ordered” states and the transition time is significantly shorter than 
the average residence time ( t r ) the system stays in the vicinity of one of the ordered states. The value 
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Figure 17: Time-dependence of the C strategy frequency in the sublattice X for the hawk-dove game on square lattice with 
nearest neighbor interactions. The fluctuations are smoothed by averaging over 100 MCS. The Monte Carlo simulation is 
performed for T = 1.5, R = 1, S = 0.5, P = 0, K = 0.5, and L = 10. 


of the average residence time can be estimated by counting the state reversals during a sufficiently long 
Monte Carlo simulation that gives — 3 x 10 4 MCS for the parameters used for results plotted in Fig. [171 
Figure [T5] shows that t av increases exponentially with the linear size for the values of payoffs and noise level 
given in the caption of Fig. 1171 Consequently, for sufficiently large sizes we may find the system staying in 
one of the macroscopic states during the whole sampling time that may exceed many years or even the age 
of the Universe. The selection of one of the ’’stationary” states depends on the initial conditions and the 
sequence of random numbers used in the simulation. Although the average residence time depends on other 
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Figure 18: MC data for the size dependence of the average reversal time obtained on a square lattice for T = 1.5, R = 1, 
S = 0.5, P = 0, and K = 0.5. 


parameters (payoffs and noise) the spontaneous symmetry breaking described above remains valid for other 
values as well as for many other systems in the limit N —> oo. This feature is contrary to the prediction of 
the Boltzmann distribution suggesting the presence of both types of ordered microscopic states with equal 
probability. In other words, the infinite system is nonergodic because the time and ensemble averages would 
give different results below the critical transition point for h = 0. An identical behavior occurs for the 
anti-ferromagnetic lattice models (for h s = 0) as well as for the Po tts m odels. This feature has rigorously 
been studied for a long time in the theory of stochastic phenomena 1871. 

It is emphasized that the above mentioned discrepancy vanishes if the equivalence of the ordered struc- 
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tures is broken by applying h =/= 0 for the ferromagnetic Ising models. In that case it is convenient to consider 
the zero-field magnetization as an average value in the limit h —> +0. A similar trick can be applied for the 
anti-ferromagnetic Ising model ( h s —>■ +0) and also for the Potts models. 

7.2. Mean-field theory 

In most of the cases the mean -field theory gives an adequate qualitative prediction about the system 
behavior. Bragg and Williams jl57l | suggested expressing the contribution of interactions via the introduction 
of an effective magnetic field h e ff characteristic to the average magnetization (m = (er x )) in translation 
invariant systems. For a d-dimensional hyper-cubic lattice 

h e ff = zJm (103) 

where z = 2d is the number of interacting neighbors. 

In the presence of a magnetic field h + h e f / the average magnetization for a single spin is determined by 
the Boltzmann distribution as 


m = tanh \{h + h e ff)/K\ = tanh \{h + 2 dJm) / K]. (104) 

This implicit equation has a trivial solution m = 0 if h = 0. This (paramagnetic) solution is the only solution 
above a critical noise level ( I\ > K c = 2 dJ). In the opposite cases (K < K c ) two additional equivalent 
(symmetric) solutions are found as illustrated on the left plot of Fig. [191 Both solutions tend towards a 




Figure 19: (Color online) Prediction of mean-field theory for the magnetization m vs. K/2dJ for the Ising model at h = 0 
(left) and h = 0.01 (right). Solid (dashed) lines denote attractor (repellor) solutions. Colored regions refer to three basins of 
attraction. 


homogeneous state ( m(K ) = ±1) and m 2 {K) vanishes linearly when K —> K c from below. More precisely, 


m{K) 


3 (K c - I<) 1/2 
Kr 


(105) 


in the close vicinity of K c . 

In the presence of an external magnetic field (h > 0) the typical solutions are illustrated in the right plot 
of Fig. [ini Notice that Eq. (llOdl) has also three solutions at low values for K while in the opposite limit the 
solution becomes unique. Contrary to the zero-field case now there is a solution that varies smoothly when 
K is increased from zero to oo. 

It is emphasized that the thermodynamically stable solution is distinguished by the extremum principle at 
the level of one site approximation. In that case the one-site probabilities are expressed as pi(+l) = (l+m)/2 
and pi(— 1) = (1 — to)/ 2 and from the resulting thermodynamical potential $(to) Eq. (I104D can be deduced 
by the extremum condition dQ(rn)/dm = 0. Direct comparison of the values of 'F(to) for the different 
solutions justifies the preference of the upper (thicker solid) curve in the right plot of Fig. [THl 
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The dynamical mean-field equation can also be used to determine the dynamical stability of the above 
solutions. The corresponding equation of motion summarizes the contributions of spin reversals for the 
Glauber type evolutionary rules, that is, 


dm 

~dt 


1 + m 1 

2 ^ g(4dJra+2/i )/K 


1 — 771 1 

2 \ _|_ g— (4:dJm-\-2h)/K ’ 


(106) 


Evidently, in the stationary state this equation becomes equivalent to Eq. (ITotI) . At the same time Eq. m 
allows us to determine the basin of attractor for each solution as denoted in Fig. [TH] by different colors. 
Namely, if the system is started from a state with a magnetization ?7i(0) then m(t) evolves vertically towards 
the corresponding attractor. For high noise levels the system always develops into the only stationary 
solution. For h = 0 and below the critical noise level ( K < K c ). however, the final state depends on 
the initial state, that is, if m(t = 0) > 0 then the system develops towards the positive m(K) stationary 
solution. For K < I\ c the trivial solution m = 0 is a separatrix. In the presence of an external field 
the separatrix is distorted as illustrated in Fig. 1191 Consequently, when varying the initial magnetization 
the dynamical mean-field theory allows the system to remain in both ferromagnetic states contrary to the 
extremum principle favoring only one showing continuous variation in m(K). 

Thermodynamically unstable (or meta-stable) states can be observed very frequently in the nature. 
For example, after a suitable variation of temperature (or other thermodynamical quantities) materials 
can evolve into (or remain in) a meta-stable phase like the supercooled liquids and gases. The disordered 
phase in alloys (cooled down fast from a high temperature) remains present for a long time despit e th e 
fact that phase segregation is favored thermodynamically at low temperatures ( e.g ., Cu-Au alloys) (163]. 
The magnetic hysteresis represents another example where the reversal of a weak magnetic field is not 
accompanied directly by the reversal of m. Similar phenomena are p rese nt in biological and social systems 
when varying the payoff parameters or dynamical rules in time 18 lam- 

In the light of the above phenomena the results of the mean-field approach can be interpreted as a 
message that the mean-field conditions do not support the system to achieve the thermodynamically stable 
state. Such a situation can occur in a social system where the players select z co-players at random from 
the whole population and determine their strategies in the spirit of the Glauber dynamics. Later we will 
show that for short-range interactions the system can evolve into the thermodynamically stable state due 
to the enhanced role of fluctuations. 

In fact, the main shortage of the mean-field theory is related to neglecting the fluctuations. This is not 
dangerous when there are a large number of neighbors (here for z = 2d 1) and the law of the large 
numbers ensures small variance in h e ff. This is the reason why mean-field type behavior is expected in 
many spatial multi-player systems if the spatial dimension d exceeds a critical value, that is, if d > d c . 

In a low-dimensional system the fluctuations affect significantly the system behavior. Due to this effect 
the prediction of mean-field theory is wrong in the one-dimensional systems (with short-range interaction) 
that do not have long-range ordered arrangements at finite temperatures. Otherwise the mean-field theory 
also gives a qualitatively good picture about the phase transitions in more complicated systems. 


7.3. Series expansions and duality 

The prediction of mean-field theory can be improved by the application of low and high noise series 
expansions allowing us to extract additional properties of the Ising model. Now we study the traditional 
ferromagnetic Ising model without magnetic field (h x = 0) on a square lattice with uniform interactions 
Jxy = J between the nearest neighbors in order to illustrate the duality that is an inherent symmetry of 
this system. 

In the low noise limit (K 0) the partition function [see Eq. (1731) ] can be approximated by summing 
the contribution of those states that are present with the highest probabilities for a finite N as 

Z = e 2NJ/K [1 + Ne~ A ' 2J/K + 2Ne~ e ' 2J/K + •••]. (107) 

Here the first term gives the contribution of the ordered ’up spins’ state (a x = +1); the second one comes from 
states where only a single spin is reversed inside the ordered arrangement, and the third term summarizes 
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the contribution of states where only two nearest-neighboring spins are flipped. The argument of the 
corresponding exponential functions reflects the interfacial energy of the island of spins reversed. 

Similar expressions can be evaluated for other lattices even in the presence of a homogeneous mag¬ 
netic field. The resulting expression can be used to derive analytical approximations for thermodynamical 
quantities ( e.g. magnetization, energy, specific heat) in the low noise/temperature limit. 

The high noise series expansion is based on the following transformation suggested by van der Waerden 


1911] , Accordingly, the contribution of a spin-pair to the Boltzmann factor is divided into two parts as 


e -z axa y = cosh (J/K)[ 1 + a x (j y tanh (J/K)} 


(108) 


where the second term goes to zero if K —> oo. Using this expression the partition function obeys the 
following form: 

Z = ^ = [cosh ( J/K)] 2N ^ n [1 + a x&y tanh (J/I\)] (109) 

<y {x,y) <y {x,y) 

where the product runs over all possible nearest neighbors. The product can be expanded out into 2 2N terms, 
where most of them [e.g. u x a y tanh ( J/K) or (<J x a y )(a x ia y i) tanh 2 (J/I\)} have vanishing contribution to 
the partition function after the sum is over the spin configurations. Exceptions are those products of spin 
pairs where each a x has an exponent of 2 or 4 when their value is one for all the 2 N spin configurations. 
The latter constellations of different pairs have a general geometrical feature, namely, the connected a x cr y 
nearest neighbors form a closed loop (or collection of loops) as illustrated in Fig. [2D] 




Figure 20: (Color online) Leading pairs of contributions to the partition function for the low and high-noise series expansions 
on the square lattice. For the low noise case (left) islands of down spins are indicated by open bullets in the sea of up spins 
(closed bullets) while thick (red) edges denote spin-pairs contributing to the interfacial energy. In the high noise limit (right) 
thick (blue) edges denote loops of neighboring spin-pairs giving non-zero contribution when summing over the spins. 


As a result, the leading contribution to the partition function can be written as 
Z = 2^ [cosh (J/K)} 2N [ 1 + TV [tanh (J/K)} 4 + 2V[tanh (J/K)} 6 + 


( 110 ) 


Comparing the expressions between the rectangular brackets in Eqs. (110711 and (11101) indicates their similarity 
that remains valid for arbitrary hig h or der of the approximati ons as justified by combinatorial methods 
detailed in the reviews by Wannier 192 1. Newell and Montroll jl59| , and Dornb |193i |. This feature serves 
as a basis for the so-called duality transformation mapping the low-noise behavior at J/K onto a high noise 
Ising system with J/K' on the square lattice if 


e ~2 J / K = tanh ( jj k >^ _ 
40 


( 111 ) 




















This formula reflects that the square lattice Ising model is self-dual, that is e~ 2J ^ K = tanli ( J/K ). I n sh ort, 
the dual of the dual is the original system. This symmetry was exploited by Kramers and Wannier [194 ] in 
the first exact evaluation of the critical temperature. 

Before going to the discussion of the critical transition we briefly mention several properties related to 
the high-noise series expansions or duality. First we emphasize that the concept of duality remains valid for 
other two-dimensional lattices. That means, for example, that the low-noise behavior on a triangular lattice 
can be mapped into the high noise behavior on the honeycomb lattice as detailed in the above mentioned 
reviews. Notice furthermore that both the high- and low-noise series expansions can be performed when the 
coupling constant depends on x therefore it involves the possibility that the low noise behavior of an Ising 
model with inhomogeneous coupling constant can be mapped onto another inhomogeneous Ising model at 
high noises (disregarding the irrelevant prefactors in Eqs. (11071) and (11101) 1 [l95l |. 

Notice furthermore, that the high-noise series expansion reflects the relevance of loops formed by inter¬ 
acting spin pair throughout the connectivity structure. So this approach indicates directly the absence of 
corrections for the loop-free networks, namely, on the one-dimensional lattice with open boundary or on the 
tree-like structures including Cayley trees. As a result, the partition function of the Ising model on loop-free 
networks can be given as 

Z = 2 n [cosh (J/K)] Nl (112) 


where A) is the number of links (connected spin pairs) in the given structure. For example, Ni = N — 1 
for the one-dimensional lattice where the vanishing correction ([2sinh J/K] N ) of the periodic boundary 
condition can also be evaluated as it comes from the only loop containing all the sites/edges. Notice that 
Ni = N — 1 for all the tree-like structure (independent of the degree distribution). 


7 -4- Critical phase transitions on lattices 

In the quantitative analysis of the critical pha se tr ansition(s) from the disordered state to the ordered 
one, several quantities play fundamental roles 196 . Il62 1 . For example, in a ferromagnetic system the average 
magnetization m is defined as 

m=—J2(<7 x ) h ^ +0 (113) 


where the sum runs over all sites x of the lattice. In general m is considered to be an order parameter where 
(• • denotes long-time average in the stationary state in the presence of an external magnetic field 

with strength h —> +0. For any finite value of h the probability of the states of disfavored magnetization 
is suppressed in the Boltzmann distribution in the limit N —> oo. Using this approach we can avoid the 
difficulties related to the nonergodicity mentioned above. 

In the equivalent evolutionary potential games the above order parameter can also be used for the 
quantitative analysis of the ordering by substituting a x = +1 ( a x = — 1) if s x = D ( s x = C). 

For the analysis of an anti-ferromagnetic system the lattice sites are divided into two disjoint sublattices 
(x £ X and y £ Y ) as discussed in Sec. [I] This anti-ferromagnetic system can be transformed into a 
ferromagnetic model by substituting J — > — J and a y —> —a y Vy £ Y. This transformation is accompanied 
by a transfer of the homogeneous magnetic field h into a staggered magnetic field h s and vice versa. In the 
absence of magnetic fields this transformation explains the equivalence of the corresponding order-disorder 
transitions. 

For these sublattice ordered spatial structures the anti-ferromagnetic order parameter is usually defined 


as 



y" (<7x)h s -y+0 - Wy)h a ^+0 

x£X y£Y 


(114) 


where the averages are evaluated in the limit h s —> +0. Henceforth we will not denote these assumptions. 

In Fig. (21]Monte Carlo data of the order parameter are indicated by symbols when varying the noise level 
(temperature) for an evolutionary social dilemma game with nearest neighbor interactions on the square 
lattice. Data obtained for Glauber dynamics at parameters T = 1.5 and S = 0.5 are identical to the case 
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Figure 21: (Color online) Monte Carlo results for the order parameter m as a function of noise for evolutionary potential games 
with Glauber dynamics on square lattice at several values of payoff parameters: in the notation of social dilemmas: T = 1.5, 
S = 0.5 (boxes); T = 1.4, S = 0.3 (diamonds); T = 0.5, S = —0.49 (open circles); and T = 0.5, S = —0.499 (closed circles). 
The exact solution for the corresponding Ising model is denoted by a solid line. Dotted (blue) and dashed (red) lines illustrate 
the prediction of the cluster variation method for the levels of one- and two-site approximations. 


of an anti-ferromagnetic Ising mod el w ithout the external magnetic fields (h = h s = 0) for which the exact 
solution was obtained by Onsager [l58| as 

1 

(115) 

where K c = 2 J/ In (1 + V2). The results illustrate clearly how the system undergoes an order-disorder 
transition at a critical noise level. As mentioned above, this value of K c can also be evaluated from the 
duality relation by substituting K = K' = K c into (11111) . For low noises (K < K c ) the order parameter 
varies from 1 monotonously and continuously to zero when approaching the critical point and remains m = 0 
if K > K c . Similar transition can be observed for other evolutionary games (where T — 1 = 5) corresponding 
to the zero-field Ising model. More precisely, the m(K/K c ) functions coincide and K c is proportional to 
|T — 1 + S\. 

It is well known that for the anti-ferromagnetic Ising model the main features of the critical transition 
remain unchanged in the presence of an external magnetic field h whereas the critical noise level decreases 
as illustrated by data obtained for T = 1.4 and S = 0.3. This behavior is related to the fact that h does not 
influence the total energy of the two equivalent ordered spin arrangements. Similar effects can be observed 
for the ferromagnetic Ising model when a staggered magnetic field (h s ) is switched on for h = 0 due to 
the above-mentioned intimate relationship between the ferromagnetic and anti-ferromagnetic Ising models. 
On the contrary, the application of h to the ferromagnetic systems smoothed the abrupt variation of the 
order parameter as illustrated by open and closed circles in Fig. [21] These Monte Carlo data show how the 
magnetization tends towards the exact solution if the corresponding magnetic field goes to zero. Notice, 
furthermore, that m remains positive for arbitrary noise levels while m —> 0 if K —> oo. In addition, the 
symmetries prescribe that m(—h) = — m(h ). Evidently, similar behavior occurs for the anti-ferromagnetic 
system in the presence of a staggered field. 

Notice the excellent coincidence between the analytical and numerical results we have obtained for 
T = 1.5 and S = 0.5 while Monte Carlo data (with m < 0.5) are missing in the close vicinity of the critical 
point. The latter deficiency is a direct consequence of the technical difficulties related to the general features 
of this type of critical transitions. Namely, the average value of the order parameter vanishes algebraically, 
more precisely, 

to oc ( K c — Kf , (116) 


l<c 


1 - [sinh (ln(l + v^)—£)] 
K 


-4 


42 










with /3 = 1/8 if K —> K c from below for all the cases (T > 1 and S' > 0) when the chessboard-like ordered 
strategy arrangement transforms into a disordered (m = 0) state as illustrated in Fig. [22] In the absence of 



K c -K 


Figure 22: Log-log plot of the order parameter as a function of K c — K for results plotted in Fig. 1211 


the exact two-dimensional solution for h ^ 0 the latter universal feature was suggested by Griffiths [1971] and 
justified by Rapaport and Domb fl98l ] who used approximation methods in the investigation of the transition 
from the disordered (paramagnetic) state into the anti-ferromagnetic structure in the presence of a magnetic 
field. The exploration of this universal behavior in the critical phase transitions was also motivated by the 
experiments. In fact, the first critical exponent (or divergency) was observed when measuring the specific 
heat as a function of temperature in the absence of magnetic field. Within the context of thermodynamics 
the different versions of the specific heat are described as the second partial derivative of the suitable 
thermodynamic potential and denoted as 


c h oc | K c - K\ c 


(117) 


where the index h refers to fixed magnetic field for the anti-ferromagnetic systems. For the ferromagnetic 
system h = 0 is required otherwise the external field h suppresses the divergence of Ch in parallel with the 
elimination of the power law behavior in the magnetization as shown in Fig. 1211 Accordingly, instead of the 
power law divergency, one can observe a peak in the A'-dependence of the specific heat at the vicinity of 
K c and the height of this peak decreases if h is increased. Evidently, the internal symmetries of the Ising 
model imply that the application of a staggered magnetic field ( h s ) results in similar consequences for the 
transitions from anti-ferromagnetic to paramagnetic phase. 

In agreement with the features mentioned above and with the fluctuation dissipation theorem [see 
with further references therein] the magnetization depends sensitively on the external magnetic field at the 
critical point (K = K r ), namely, 

m(x\h\ 1/s . (118) 

In the understanding of the universal features of the critical behaviors, the analysis of the correlation 
functions played crucial roles. The general versions of correlation functions are defined in the translation 
invariant stationary states as 

3 (, ’ j) (x,t) = (n^(t')n { x J l y (t +1')) - (nW(r))(n^ y (t + 1')) (119) 

where i and j denote strategy labels ( i,j = 1,... ,n), (• • •) refers to averaging over all sites y and time t'. 
Notice that here we use the vector notation of sites in the argument of the correlation functions although in 
many numerical analyses x refers to horizontal (or vertical) spatial distances, n^ is the extended definition 
of the occupation variable at time t and site x, namely, 




1 for s x = i , 
0 otherwise 


( 120 ) 
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For sublattice ordered structures it is convenient to prescribe that x £ X and y £ Y. For short range 
interactions the one-time correlation function decreases exponentially with the horizontal or vertical distance 
as 

</*’ J )(x, 0) ~ e"^ 1 (121) 

if |x| —> oo. The one-site correlation function shows similar behavior, namely, 

g ( ~ i ' j \0,t)~e T v ( 122 ) 


if t 1. Both the correlation time Ty and the correlation length j depend on the noise level K as well 
as on the potential parameters. In the Ising model for h = 0 we can omit the indices (i and j) and the 
corresponding quantities diverge as 

£ oc \(K C — K)\~ v± (123) 

and 

roc \{K C -K)\~ V \\ (124) 


with v± = i/|| = 1 when approaching the critical point. We note that characterizes the so-called critical 
slowing down in the time-dependence of correlations. Beside it, the intensity of the fluctuation of order 
parameter given as 

r -i 2 


Xm=N 


( 


(Y a ^-Y' 


(125) 


also diverges. More precisely, 


Xm oc \(K C - K)\ 7 . 


(126) 


In the stationary state at the critical point the correlation function decreases algebraically, that is, 


#(x,0) oc |x| 2 d ri± . 


(127) 


All these features jointly cause serious technical difficulties in Monte Carlo simulations if one wishes to 
determine the order parameter with a reliable accuracy in the close vicinity of the critical point, because 
both the system size and sampling time should be chosen to be significantly larger than the correlation 
length and time. 

In fact the exponents introduced above are not independent. One of the main results of the statistical 
physics in the 20th century is the explanation of the universal behaviors in the critical phase tr ansit ions. 
Using t he s caling hypothesis and the renormalization-group techniques, introduced by Kadanoff [200l | and 
Wilson [201 ], relations have been explored between the above defined critical exponents. One form of this 
idea assumes that the singular part of the thermodynamic potential <f> (K, h) in the close vicinity of the 
critical point K c is a generalized homogeneous function of s = \K — K c \/K c and h, e.g., 


<S> s (\ ah h,\ a *e) = \$ s (h,e) (128) 

as surveyed briefly by Stanley |202| . This concept implies that two parameters (ah and a e ) determine the 
values of critical exponents. The validity of the scaling hypothesis is confirmed by numerous experiments 
and theoretical calculations by observing data collapse when a thermodynamic quantity of two variables 
(e.g. m(h,K — K c )) forms a single curve, if we use suitable scales in a two-dimensional plot. The latter 
feature reflects a universal behavior and is valid in the close vicinity of the critical point for the systems 
belonging to the Ising universality class. 

The concepts of renormalization group techniques helped us understand the phenomena and relevant 
conditions resulting in the universal behaviors (for a brief survey of the essence and history of this approach 
see the review by Fischer |203| l. This approach distinguishes relevant (e.g., spatial dimension d, number 
n of possible states, and additional symmetries) and irrelevant (lattice structure) quantities and gives an 
adequate description and explanation of the universal behavior at the critical point (and also in its close 
vicinity). These investigations have justified that the diverging quantities tend asymptotically towards a 
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exponent 

d =2 

d =3 

d=4(MF) 

a 

0(log) 

0.110(1) 

0 

p 

1/8 

0.34265(3) 

1/2 

7 

7/4 

1.2372(5) 

1 

<5 

15 

4.789(2) 

3 


1/4 

0.0364(5) 

0 

v± 

1 

0.6301(4) 

1/2 


Table 1: Critical exponents of the d-dimensional Ising model. 


power law behavior characterized by the same exponent on both sides of the critical point. For example, 
when approaching the critical point £ ~ C\(K c — K) 1 if K < K c and £ ~ — X c ) 7 if K > K c (C\ ^ C 2 ). 

This is the reason why we used the same notation for the exponents below and above the transition point. 

The renormalization group techniques have confirmed quantitatively the universal scaling hypothesis 
between an upper and lower spatial dimension d. The predictions of the scaling hypothesis are simple 
consequences of the properties (11281) preserved by the Legendre transform of &(K, h) and its derivatives. 
Since the critical exponents are directly related to ah and a E , one can derive scaling laws by eliminating 
these variables. Using this method three scaling laws can be derived: 

a + 2/3 + 7 = 2, (129) 

a + P(6 + l) = 2 , (130) 

(2-r]±)v ± = 7 . (131) 


which are independent of the spatial dimension d, whereas the fourth, so-called hyper scaling law , indicates 
directly the relevance of the spatial dimension as 

2-a = dv x . (132) 


Due to these relations only two of the six static exponents are independent. For the practical identification 
of this class of universal behavior in any models one needs to check the values of only two static exponents. 

For the systems belonging to the Ising universality class, the values of most relevant exponents as a 
function of the spatial dimension d are summarized in Table [TJ The reader can fin d a c omprehensive list 
of the th eoret ical and experimental results in the review by Pelissetto and Vicari |l68l | and in the book 
by Odor [2041 ]. The one-dimensional results are missing here due to the absence of the critical transition. 
For d > 4 the critical behavior and exponents are well described by mean-field theory. Furthermore, the 
two-dimensional results are extracted from the exact results and the value of a = 0 refers to logarithmic 
divergence (Ch oc — logs) in the specific heat. 

Finally we emphasize that a similar universal critical behavior can even occur in other non-equilibrium 
models where the probability distribution of the microscopic states differs from the Boltzmann distribution. 
For example, Perez et al. [203 ] have reported this type of universal critical transition in an Ising model, 
where the spins are reversed simultaneously as in the stochastic cellular automata surveyed by Wolfram 

f20i . 


The critical phenomena have also been explored for the Potts models on different lattices and graphs. 
Most of our knowledge comes from a wide scale of analytical (approximation) methods or numerical tech¬ 
niques and are justified by experiments. As it is known, the mean-field descr iption gives a qualitatively 
correct picture of the phase transition in the Ising model. Kihara et al. [l80 ] found that the mean-field 
approach predicts first-order phase transition for all n > 2. Subsequently, the more sophisticated methods 
have clarified that the order-disorder transitions are of first-order for n>4atd = 2;n>3atd = 3; and 
n > 2 at d = 4 in homogeneous spatial systems 


1821 - 


Contrary to the mean-field prediction, the order-disorder transitions are continuous for all the two- 
dimensional lattices if n = 2, 3, or 4. These critical transitions exhibit power law behavior in several 
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exponent 

n = 2 

n = 3 

n = 4 

a 

O(log) 

1/3 

2/3 

p 

1/8 

1/9 

1/12 

7 

7/4 

13/9 

7/6 

<5 

15 

14 

15 

V 

1/4 

4/15 

1/2 


1 

5/6 

2/3 


Table 2: Critical exponents vs. n for the two-dimensional Potts models. 


quantities in the close vicinity of the transition point as detailed above. The corresponding critical exponents 
are listed in Table [2] where the values in column n = 2 are equivalent to the critical exponents of Ising model 
given also in Table [L] Evidently, the listed values of the critical exponents satisfy the relations expressed by 
unsD-ra as these are derived from general scaling laws. 

It is emphasized that the above robust behavior occurs in homogeneous spatial systems where n individual 
strategies or degenerated ground states can be distinguished. Deviations from these universal behaviors 
appear when the spatial systems become inhomogeneous or the connections between the players are described 
by different graphs [207j| . 

7. 5. Ordering in other relatives of the Ising model 

The decomposition of payoff matrices for n strategies has highlighted the relevance of those basis games, 
denoted as d(p, q), that possess Ising type interactions between the pth and gth strategies. These models 
have not yet been investigated systematically. The Monte Carlo simulations on a square lattice with logit 
dynamics indicate an Ising type ordering process when increasing the noise level K as it is illustrated in 
Fig. [23] for three values of n. For all the three interactions d(l,2) the system prefers the homogeneous 
s x = 1 or s x = 2 state in the low noise limit and (p± — p 2 ) —» 0 if K approaches K c from below. The 
preferences of the strategies 1 and 2 can be observed for K > K c because of the formation of their cluster 
when pi = p 2 > p n (for n > 2). However, in the limit K —> oo the fluctuations suppress these correlations 
and all the n strategies are present with the same probability (1 /n). 



Figure 23: (Color online) Strategy frequencies as a function of K for spatial evolutionary games if the interaction is defined by 
d(l,2) [see Eq. (1411) 1. The Monte Carlo data of the 1st (2nd) strategies are indicated by closed (open) symbols for n = 2 (red 
o) ? 71 — 3 (A), and n = 4 (blue boxes). The frequencies of the third and fourth strategies for (n = 4) are indicated by (blue) 
pluses and crosses while the symbol v shows the frequency of the third strategy for n = 3. 


The first numerical investigations 31, 59), SflB] support the theoretical expectations predicting Ising type 


critical transitions at K c (n) for n = 3 and 4. The preliminary Monte Carlo results show that K c [n) decreases 
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if n is increased and the phase transition becomes a first order one if n exceeds a threshold value. 

In the above-mentioned systems the critical phase transitions are smoothed if the interaction is perturbed 
by switching on a self-dependent component or an additional coordination type interaction that can prefer 
the first (or the second) strategy. The Pi(K) functions in Fig. [24]are obtained for an interaction that results 
in the dominance of strategy 1 in the low noise limit. In this plot we were able to illustrate the Monte Carlo 
data by lines due to the absence of the fluctuation enhancement in the region of smooth transition. 



Figure 24: Strategy frequencies vs. noise level K if the interaction is given by A = d(l, 2) + 0.1 ■ d(l, 3) in the three-strategy 
multi-agent evolutionary potential game. 

As mentioned previously, the symmetric potential matrix can possess two equivalent preferred Nash 
equilibria [(p, q) or (q,p) if max(I4,) = V pq = V qp \ even for n > 2. On the square lattice these interactions 
yield two equivalent sublattice ordered strategy arrangements when the players select strategies p and q 
in the sublattices X and Y or conversely if A' —> 0. When increasing the noise level K these systems 
exhibit an order-disorder phase transition at K c . The strategy frequencies in the sublattices show a similar 
A'-dependence plotted in Fig. [231 Contrary to the case of equivalent homogeneous ordered states, the equiv¬ 
alence between the sublattice ordered strategy distribution is not destroyed if the interaction is perturbed 
weakly because their effects are identical on both sublattices. Consequently, the Ising type sublattice order¬ 
ing takes place in a wide region of the payoff parameters. The latter phenomenon is resembling the effect 
of the external magnetic field h on the anti-ferromagnetic ordering if its value does not exceed a threshold 
value (see Figs. flGl and ITU. 


7.6. Critical phase transitions on networks 

The investigation of the Ising model on different networks also has a long history as detailed in the 
review by Dorogovtsev et al. [207l |. First we discuss the simplest cases where the players are distributed on 
a Cayley tree or Bethe lattice that is considered as an infinite Cayley tree without its periphery. Despite the 
strong relationship between these connectivity structures, the Ising model exhibits fundamentally different 
behaviors on them a s it w as emphasized in early works of Eggarter 111], Muller-Hartmann and Zittartz 
209 . Wang and Wu [2JUJ. 

The first important difference is due to the fact that a relevant portion of the nodes of the Cayley tree 
belongs to the periphery (the leaf nodes in graph theory terminology) where the players have only one 
neighbor and the presence of this boundary affects the macroscopic behavior significantly, due to the large 
num b er of the se sites. On the other hand, this feature is utilized in the determination of the exact solution 
[llll lll2L 1211 |. The absence of the ordered states on the Cayley tree for any finite values of K is related 
to the fact that when removing a single edge this structure is divided into two independent parts where 
opposite ordered structure can be formed. Due to this feature, arbitrary deviation in the magnetization can 
be achieved by generating a single point defect that is always present at finite K in the sufficiently large 
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systems. A similar reason explains the absence of ordered strategy arrangements (and phase transition) in 
the one-dimensional lattice where one_of the long-range ordered phases occurs only for K = 0. 

Miiller-Hartmann and Zittartz [209| have shown that the thermodynamic potential (free energy) on the 
Cayley tree becomes a nonanalytic function ofthe magnetic field h below a critical temperature Kbp given 
by the pair approximation applied by Bethe |l40| and Peierls jl55| . More precisely, the leading part of the 
nonanalytic behavior is proportional to h K where the exponent k increases monotonously from 1 to oo as the 
temperature goes from 0 to Kbp■ This phenomenon is accompanied by a diverging magnetic susceptibility 
below the transition point. Melin et al. [212| have shown that this unusual behavior may be related to the 
large number of metastable states that are stable with respect to single-spin flips. 

The use of Bethe lattice eliminates the boundary layer and assumes translation invariance ( i.e . the 
equivalence of sites). Using this assumption we can derive analytical results that predict mean-field type 
order-disorder phase transition if K is increased. Due to the absence of loops in this structure the two-site 
cluster variation method (called also Bethe-Peierls or pair approximation) reproduces the exact result as 
indicated in Fig. (23 where the MC data are obtained on a random regular graph for large si zes as the MC 
simulations cannot be performed on an infinitely large Bethe lattice. In network analysis 213, l214j and 



Figure 25: Order parameter m as a function of noise K for evolutionary hawk-dove game with Glauber dynamics at payoff 
parameters T = 1.5, S = 0.5 if the players are located on a random regular graph of degree 4. Xs indicate MC data for 
N = 2 X 10 6 . The solid line denotes the prediction of the two-site cluster variation method for the Bethe lattice and the 
Onsager results (on square lattice) are illustrated by dotted line. 


graph theory [55| it is well known that the average length of the shortest loops for the random regular 
graphs increases with In AT. Consequently, this approach can give an adequate description of the system 
behavior for those cases where the sparse and long loops of the connectivity structure do not modify the 
behavior relevantly. This approach is capable of explaining more striking differences in those phenomena 
where a random regular graph is substituted for the square lattice in a three-strategy evolutionary game 
f215l . 


In the last decades the statistical analysis of different networks initiated the reinvestigation of traditional 
lattice models on a wide scale of complex networks. Most of these networks are claimed to provide a better 
description of the connections that have emerged in social interactions, neural systems, communication, 
transportation, etc. Some of these network models are capable of describing a continuous transition from 
a spatial lattice to random or random regular graphs, while others create novel features not involved in 
the traditional concept of lattices. Now we briefly outline what happens to the Ising and Potts models on 
several types of networks. 

The so-called small-world effect takes place in networks where the average shortest distance between two 
nodes is proportional to In TV. The Bethe lattices, Cayley trees, and random regular graphs possess this 
feature. Watts and Strogatz [21fil | proposed a way of building the small-world effect into spatial lattices. In 
the first model they considered a one-dimensional chain (ring) with first- and second-neighbor interactions. 
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exponent 

v > 5 

v = 5 

3 < v < 5 

a 

1st order 

0(log) 

(5 - v)/(v - 3) 

p 

1/2 

1/2* 

l/(v-3) 

7 

1 

1 

1 

<5 

3 

3* 

v-2 


Table 3: Critical exponents vs. v for the Ising models on scale-free graphs. * refers to logarithmic corrections described in 



The small world feature is achieved by removing qN connections from the lattice and connecting one of their 
end points to another site chosen at random (0 < q <C 1). This method does not modify the average number 
of neighbors though it destroys the regularity. At the same time, it can be applied for any d-dimensional 
latti ces. The number of neighbors is conserved at each site by the method proposed by Newman and Watts 
[217i | where connections are interchanged between pairs of connected sites selected randomly. Soon after the 
network models had been published, different phenomena were investigated on the resultant structures. 

Studying the one-dimensional Ising model with the above rewired la ttice s or with a lattice obtained by 
adding a portion q of new links, Barrat and Weigt |218| and Gitterman j219t | have observed the appearance 
of ordered structure at finite noise levels ( K < K c (q)) in sufficiently large systems. It is found that K c 
increases monotonously with q and K c oc — z/\og{q) for low values of q where z (z > 4) is the number of 
neighbors in the starting one-dimensional chain. They also reported a mean-field type phase transition at 
K c . 

These types of small-world structures were reinvestigated by Herrero [22l| and Chatterjee and Sen 221 1 
when the initial structure was a square (d = 2) or cubic (d = 3) lattice. Their analyses have confirmed the 
presence of the ordered state at low temperatures in agreement with the expectation. The increase of K c 
was similar to those found for the one-dimensional cases, that is, SK C oc l/log(<?) and the phase transitions 
were of mean-field type. 

Up to now we have studied the Ising models on regular or quasi-homogeneous networks where the 
degrees of nodes were close to their average value. Fundamentally different behaviors were reported by 
Dorogovtsev et al. [222| and Leone et al. ]223j | who studied the Ising model on random graphs with a degree 
distribution having a fat tail. Most of these investigations are performed on networks having a power-law 
degree distribution p(z) oc z~ v for large number (z) of neighbors. Using different approaches these authors 
studied the general features of the ferromagnetic phase transitions for different values of v. The analyses 
showed mean-field type phase transitions when v > 5. These approaches indicated non-trivial critical 
exponents for 3 < v < 5 as listed in Table [3] Note that behavior of the fluctuation is not affected by the 
value of v within this region, while a, /?, and 7 depend on v in a way that breaks the scaling relations (11291) 
and (11301) . 

For v < 3 the role of the high degree nodes (hubs) becomes relevant in the limit N —> 00 and the system 
stays in a ferromagnetic state for arbitrary noise level K. Within this range of v. however, a size effect 
can be observed because the value of N limits the maximum of z, too. As a result, for most of these finite 
systems the Ising model exhibits an order-disorder transition at K c oc In (N) according to t he th eoretical 


investigations mentioned above, in agreement with the results of Monte Carlo simulations [224]. Many 


realistic networks and network models belong to this class, therefore the latter feature governs the behavior 
in both biological and social systems. 

In the above small-world networks the additional random links can be interpreted as a way of introducing 
long-range interactions. Similarly to the Levy flights, the strength of long-range interactions can be weakened 
by creating the additional links with a probability p{l) oc Z~ A where l is the Euclidean distance between 
the sites to be connected and A > 0. The behavior of the resulting system is similar to those obtained by 
introducing long-range interactions, as discussed in the following section. 

2251 2261 ] who studied the Ising and 


Fundamentally different behaviors were reported by Gefen et al. 


Potts models on several fractal lattices discussed by Mandelbrot [227 ]. In the 80s the analysis was extended 


49 


















by Bhanot et al. 228 1. d’Auriac and Rammal [229}, Bonnier et al. [Hi Monceau et al. }23l[ . Some of 
the fractal lattices were generated by an algorithm producing a Sierpinski carpet. For example, a large 
hypercube of a d-dimensional lattice with n kd sites is divided into h d equal hypercubes and ( fi d — fh) of 
them are removed. The same process is repeated for the rest of the smaller hypercubes. After the fcth 
segmentation steps we get a fractal lattice characterized by a fractal dimension df = In (m) / In (n) that can 
be tuned from 0 to d. The resultant self-similar structures made these types of models attractive for the 
application of renormalization group techniques and other exact methods. A series of investigations have 
clarified that the behavior of the Ising model on these lattices depends not only on the fractal dimensions, 
but i s als o affected significantly by other topological features like the order of ramification and lacunarity 

2251 1235. 


The order of ramification 7 Z measures the number of links to be cut in order to isolate an arbitrarily large 
portion of the network. If 7Z is finite then the ferromagnetic ordering is missing at finite noise levels. Notice 
that this criterion is the generalized version of those we used above for explaining the absence of ordering 
in the one-dimensional lattice and Cayley trees. Among the fractal lattices the quasi-linear structures 
(generated similarly to the Koch-curve) [233} and some versions of the Sierpinski gasket 234 ] represent 


networks on which the Ising model has no long-range order (and phase transition) at finite noise level K. 

On the contrary, phase transition occurs if 7Z = oo in the limit N — > oo. The first papers [234} indicated 
clearly that some critical exponents depend also on the lacunarity A (235} . which quantifies the deviation 
from the translation invariance. For the above mentioned Sierpinski carpets the lacunarity A measures how 
the removed small boxes are distributed. At high lacunarity the removed boxes form a large hole, while at 
low lacunarities these boxes are distributed ” homogeneously”. 

The critical exponents of the result ant p hase transition are studied quantitatively only on a few fractal 
lattices. For example, Monceau et al. [2311 ] and Carmona et al. [ 232} studied the critical behavior on the 
above mentioned Sierpinski carpet for d = 2, h = 3, m = 8 where df = In (8)/ In (3) = 1.8927 .... Despite the 
technical difficulties their results are consistent with the suggestion of Wu and Hu [236} stating the existence 
of a weaker universality class. Accordingly, the static exponents may vary with the geometrical parameters 
of the fractal ( e.g ., df and lacunarity) and the relations between the critical exponents (I129D - (I132I) are valid 
if df is substituted for d in (11321) . In a subsequent paper Monceau and Hsiao [237 ] studied the critical 
behavior on fractals sharing the same fractal dimension but having different lacunarities. It was found that 
the long-range order at the critical point decays faster when the lacunarity is increased for a given df. 

At the end of this section it is worth mentioning that most of the above discussions are based dominantly 
on a pproximate results supported by numerical sim ulations. Very recently, however, Dembo and Montanari 
[238} . Montanari et al. [239} . and Dembo et al. [240 ] have performed a more rigorous mathematical analysis 
of the Ising and Potts models on locally tree-like structures and their results support the picture sketched 
above. 


7.7. Long-range interactions on lattices 

In ionic compounds the Coulomb interaction plays a crucial role in the formation of the microscopic 
arrange ment s of ions even in the cases when screening complicates the analysis in the presence of opposite 
charges [165] . For the ferromagnetic materials the dipole-dipole interaction is responsible for the emergence 
of a prop er m agnetic domain structure that results in an almost zero remanence magnetization in the soft 
magnets [l53} . In solid solutions the interaction between two large interstitial atoms is mediated by the 
lattice distortion and it can be well approximated by a suitable long-range interaction. This interaction 
drives the segregation process. In social and biological systems the diffusion of opinion or chemical products 
can also mediate a similar long-range interaction. For the investigation of these type of interactions, the 
formalism of the Ising model provides a good mathematical background, although huge technical difficulties 
arise from the large number of interactions to be considered in the real phenomena. In short, the analysis 
of the Ising and Potts models with long-range interactions on lattices is recently at a beginning stage. 

There are, however, several results that are worthy of a brief outline. In many cases the coupling constant 
J xy between the sites x and y is expressed by an algebraic function J xy = l/|x — y| d+A where |x — y| is 
the Euclidean distance on the d-dimensional lattice and A quantifies the type of long-range interactions. 
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Assuming a ferromagnetic interaction, the total effect (sum) of the coupling constants J xy over the whole 
system helps preserve the ferromagnetic state. This quantity as well as the total energy, however, diverges 
if A < 0 as 

1 r°° 1 

-dr (133) 
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where a(d)r d ~ 1 refers to the surface of a d-dimensional spher e of radius r. Consequently, the ferromagnetic 
state is stabilized if A < 0. Since the early work of Ruelle |241] the systematic analysis of the Ising and 
Potts models has clarified what happens for A > 0. 

Ruelle J241] proved rigorously the absence of long-range order in the one-dimensional Ising ferromagnetic 
model if A > 1. On the other hand, Dyson [242 ^ proved the existence of a phase transition in t hese systems 
if 0 < A < 1. For A = 1 the existence of phase transition was justified later by Dyson )243j | and by 


Imbrie and Newman 244 ] who described an intermediate phase where the two-point correlation function 
exhibits powe r law decay with a .A-dependent exponent below K c . Using the renormalization group approach 
Fisher et al. 2451 ] have shown that the one-dinrensional Ising model exhibits mean-field-like behavior for 
0 < A < 0.5 while non-trivial (A-dependent) exponents are found for 0.5 < A < 1 . For higher dimensions 
the analyses were focused on two regions of A where mean-field or A-dependent exponents are predicted by 
the renormalization group techniques. 


7.8. Sublattice ordered structures on lattices 

If the range of interactions exceeds the nearest neighbors on a lattice then we may face a wide variety of 
sublattice ordered states in the low noise limit. The extensive investigations of these systems were unavoid¬ 
able in solid state physics where the effects of the second- and third-neighbor interactions with different 
strengths are studied in many materials. For example, the adequate description of the oxygen ordering in 
the Cu-0 layer of the super-conducting YE^CuaCU-i compounds required to take into consideration four 
terms of interactions |246j . 

The introduction of the additional second- and third-neighbor interactions (denoted as J 2 and J 3 ) does 
not cause qualitatively different behaviors if all the coupling constants favor the ferromagnetic structure, that 
is, if Ji, J 2 , J 3 > 0. The presence of the additional terms increases the value of K c . Fundamentally different 
consequences occur, however, if the additional interactions do not support the formation of the ordered 
arrangement dictated by the first-neighbor interactions. In the latter cases numerous types of ordered 
structure may emerge. For the illust ratio n of these phenomena we study now a simple model in vesti gated 
in detail first by Binder and Landau 2471 and recently by Yin and Landau [2481 ] and de Queiroz [2491 ] . 

Consider an Ising model on the square lattice with equivalent anti-ferromagnetic interactions (Ji = 
J 2 = — 1) between the first and second-neighbors for the presence of a magnetic field h. One can assume 
that all the possible ordered structures can be built up by tiling, that is, by repeating one of the 2 x 2 
patterns shown in Fig. 1261 The resulting long-range ordered arrangements can be well described by a four- 
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Figure 26: All the possible spin arrangements within a 2 x 2 cluster of sites that are used for generating four-sublattice ordered 
arrangements on a square lattice. 


sublattice description defining the average magnetization for each sublattice. Notice that patterns within 
the second, third and fifth column of Fig. [20] can be transformed into each other by a rotation of 90° and the 
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corresponding ordered structures have the same potential (or energy) in the macroscopic system. All these 
three types of ordered structures are fourfold degenerated and involve the possibility that the corresponding 
order-disorder transitions belong to the Potts model universality class for n = 4 and d = 2. 

Contrary to the naive expectations, simulations have indicated a more complex behavior both in the 
spatial patterns and in the phase transitions. FigurcITTlillustrates four snapshots when increasing h. All these 
simulations are performed at a low noise level in order to keep the strip like anti-ferromagnetic arrangement 
(see the bottom right plot in Fig. [27]) almost free of point defe cts. The applied magnetic fields are selected 
according to the phase diagram reported by Yin and Landau 248} and represent different phases. Besides 
these four phases there exists a disordered (paramagnetic) spin arrangements stabilized for sufficiently high 
noise level (K > K c [h )) or high fields \h\ > 8. 




Figure 27: Ordered and partially ordered spin arrangements in the kinetic Ising model on the square lattice when the first- and 
second-neighbor anti-ferromagnetic coupling constants are equivalent and K = 0.3 • \J\ for different external magnetic fields 
( h ) indicated by figures inside the (50 x 50) snapshots. 


The fourfold degenerated strip-line anti-ferromagnetic structure undergoes a phase transition to the 
paramagnetic phase at K c {h = 0) = 2.08. This curious phase transition will be discussed later. In the zero 
noise limit this ordered structure is stable as long as \h\ <4. For 4 < \h\ < 6 a partially ordered structure, 
called row-shifted (2 x 2) phase, can be observed in the low noise limit. The corresponding spin arrangements 
can be constructed from the completely ordered pattern built up from one of the (2 x 2) blocks of the second 
columns in Fig. [TS] that consist of ferromagnetic and anti-ferromagnetic rows (columns) alternately. For the 
present interactions the total energy (or potential) is not changed if the anti-ferromagnetic rows are shifted 
horizontally at random. In the resultant partially ordered phase the entropy remains finite in the limit 
K —> 0. Due to the symmetries the equivalent column-shifted states can be constructed in a similar way. 
Thus within a region of h dependent on the temperature this system has two equivalent sets of partially 
ordered structures. Within this region one can observe a spontaneous symmetry breaking that is analogous 
to the formation of ferromagnetic or anti-ferromagnetic structures if the system is started from an ordered 
structure at a low noise level. This process is driven by the increase of entropy while the coexistence of the 
domains of row- and column-shifted phases is prevented by the positive interfacial energy. 

Betwee n th e strip-line and row-shifted structures, it has been possible to distinguish two additional 
phases by [248} in the phase diagram at low noises. The upper-right snapshot of Fig. [TB] illustrates the spin 
arrangement in the paramagnetic phase with short-range correlations resembling all the above mentioned 
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ordered phases. The appearance of the disordered (paramagnetic) phase is common at low noises in the 
phase diagrams between basically different ordered phases. In the present model, however, another type 
of structure (see the left-bottom snapshot in Fig. eh is observed. The latter phase can be considered as 
a poly-domain pattern of the strip-line phases where the p rese nce of interfaces is stabilized as it occurs in 
water-oil emulsions or other self-assembling systems 25(1 2511. Remarkably, the distinction of this phase 
required sophisticated techniques. 

for the 


Recently different methods have been applied by Yin and Landau [248j and de Queiroz [24S 
classification of the phase transitions in the present model. According to these investigations, the phase 
transitions from the ordered or partially ordered structure into the disordered one are continuous but not 
universal. Both transitions exhibit power law behavior with exponents dependent on the magnetic field. In 
other words, these transitions do not belong to the above mentioned universality classes of the Ising and 
Potts models. 

The above features are not unique. Similar phenomena can be obse rved for the present model when 
vary ing the ratio J2/J1 or when additional interactions are switched on [2521 ]. In the paper by Dublenych 
[253l | the reader can find dozens of patterns that can be used for building up long-range ordered structures of 
the Ising model on triangular lattices, if first- and second-neighbor interactions are taken into consideration. 
Evidently, analogous phenomena are expected in lattices for higher dimensions if the range of interactions 
is increased and also in evolutionary potential games when the number of strategies exceeds 2. 


7. 9. Frustration 

In the previous section we have faced situations when the ground state was infinitely degenerated (see 
the upper-left snapshot of Fig. [26]). Such a situation becomes natural on lattices where three spins (players) 
interact with each other via equivalent anti-ferromagnetic pair interactions (anti-coordination games) as 
denoted on the left panel of Fig. EH 254j. This phenomenon was first described by Wannier |255| who studied 
the anti-ferromagnetic Ising model on a triangular lattice. In fact, for the anti-ferromagnetic interactions, 
the basic reason is related to the presence of loops with an odd number of edges. Such topological situation 
occurs on many other lattices with a finite clustering coefficient ( e.g ., square lattice with the first- and 
second-neighbor interactions), in some curious lattice structures [e.g., the Cairo pentagonal lattice studied 
by Rojas et al. [256| ] that includes five-edge loops. 

On the contrary in social networks the length of loops may be either odd or even, consequently, the 
frustration is unavoidable on these structures for the anti-coordination games. 



© 


Figure 28: Two typical topological conditions yielding frustration. The left panel shows three players with two strategies 
(closed and empty bullets) located on a triangle for pairwise anti-coordination games denoted by double lines. The question 
mark refers to the equivalence between the two strategies for the third player. The right panel shows four players on a square 
network who play pairwise coordination (simple edge) or anti-coordination (double edge) game with each other as denoted. 
For the given three-player strategy profile the fourth player has two equivalent strategies. 


The other source of frustration is related those models where the coupling constants are chosen to be 
+J or — J at random as it is illustrated in the right panel of Fig. EH Some consequences of the frustration 
are well illustrated in Fig. EH where one can observe six different domains representing the chessboard and 
anti-chessboard arrangements of two of three strategies o n the square lattice for games equivalent to the 
repulsive three-state Potts model. Berker and Kadanoff |257l | have shown that the typical size of these 
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Figure 29: (Color online) Typical snapshot for the repulsive three-state Potts model on the square lattice at low noise levels. 


domains (or the correlation length) increases algebraically when the noise level K is decreased. Similar 
features were reported previously by Wannier [2551 ] for the Ising model on a triangular lattice. 

A direct consequence of the frustration can be the high level of degeneracy of states having maximal 
potential. If the number of these degenerate states increases with N exponentially then the specific entropy 
remains finite in the zero noise limit, contrary to the third law of thermodynamics. This discrepancy 
vanishes if the games are not equivalent, as it occurs in social systems. On the contrary, if we study a 
system with many edges of small potential differences in the dynamical graph, then the systems may reach 
the thermodynamically stable state faster. 


7.10. Effects of randomness 

For many multi-agent systems we cannot assume equivalent games representing the interactions among 
the participants. In Sec. 16.11 it is shown that if the interactions are limited to symmetric 2x2 games then 
the potential of this multi-agent system can be mapped onto an Ising model where the values of J xy and 
h x vary within a wide range. The effects of these types of randomness have been studied for four decades 
and the results are documented in the literature of spin glasses. The experimental investigations of these 
syste ms were first motivated by the unusual magnetic behavior of some metal alloys ( e.g. AuFe and CuMn) 
[258J. From a theoretical point of view these phenomena demanded the development of new concepts and 
approaches. Comprehensive de scrip tions of the models, methods, phenomena, and perspectives are given 
in the books by Mezard et al. [259l | and by Stein and Newman 12601 ]. For most of the investigations the 
effects of randomness (occurring in J xy and h x ) are considered separately. Exceptions are represented by 
the one-dimensional random Ising models 126 it . Recent trends in the research of evolutionary game theory 
offer further possibilities for the introduction of randomness via the consideration of personal features and 
the wide scale of connectiv ity st ructures. Here w e sur vey onl y a fe w phenomena t hat seem to be important 
for information processing | l262 |. decision theory 263 1 , social |260| and biological 264 1 systems. 

Spin glasses are systems combining quenched randomness and frustration. For the Ising type spin glasses 
the parameters J xy and h x are random variables and are characterized by some probability distribution 
functions. The accurate mathematical form of the probability distribution of the J xy values is not important 
in general. For technical reasons the Gaussian distribution is chosen with a mean value of zero (denoted as 
(( Jx V ))j = 0 and {{J xy Y)j_= 1 where ((• • •)) j refers to averaging over all pairs). An alternative approach was 
suggested by Domany 2651 ] for the so-called ±J model where J xy = J (J xy = — J) is chosen at random with 
a probability p (1 — p). Here we remind the reader that (( J xy ))j > 0 does not guarantee the ferromagnetic 
order in the low noise limit in the absence of magnetic field ( h x = 0). When investigating the ±J model on 
a square lattice Domany |265| has explained the existence of the ferromagnetic order in the low noise limit 
if p > p c ( p c cs 0.83 on the square lattice). It was found that the percolation of unfrustrated squares was 
responsible for the ferromagnetic state and the critical point ( p c ) is associated with the adequate percolation 
threshold value. 


54 












The spin glass model was introduced and studied by Edwards and Anderson 266] for a uniform Gaussian 
distribution ((( J X y))j = 0 and {(J xy ))j = 1). In agreement with the experimental results this model has 
indicated the presence of a ’’spin glass” phase at low noise levels ( I\ < K sg ). The main characteristic of 
the spin glass phase is the existence of a large number of free-energy valleys that are inaccessible from each 
other in the limit TV —> oo. Within the spin-glass phase the average magnetization is zero and there are 
ferromagnetic domains with different sizes and orientations. During the evolution the system remains in the 
vicinity of one of the given free-energy valleys similarly to the ferromagnetic phase where only two (ordered) 
microscopic states are d istin guished. For the quantification of the remanence of these spin arrangements 
Edwards and Anderson [266] have introduced an order parameter which is similar to the one-site two-time 
correlation function characterizing the coincidence of the ’’average” spin directions {(Jx 1 ) and (ai 2 ' ) ) in the 
vicinity of times t\ and £2 which are far from each other. If the system remains in the same free-energy valley 
then (cri 1 ^) = (c4 2 ^} = (a x ) and the spin glass phase can be characterized by the so-called Edwards-Anderson 
order parameter: 

4E(iw’» f (i34) 


qEA 


At high temperatures qEA = 0 because each spin can flip into the opposite state within a short time. 
qEA = 1 refers to frozen patterns at K = 0. Deviations from qEA = 1 arise from thermal fluctuations and 
frustration. 

A more detailed insight into the nature of spin glasses is provided by considering an exactly solvable model 
introduced by Sherrington and Kirkpatrick 12671 . It is an Ising model on a complete network (interactions 
exist between all spin pairs), therefore this model could be well investigated by the mean-field theory using 
the replica trick. For the analytical analysis of this model, [267] [268[ have justified the phase transition at 
K c = 1 from the paramagnetic state to the spin glass phase when decreasing the temperature and explained 
the magnetic behavior observed in experiments. 

In order to study the inherent relationship between two glassy states Parisi [13 has suggested using a 
more complicated order parameter: 

q af3 = J-f ^(CT*)a( Ox)p- (135) 


that quantifies the average overlapping between two frozen spin patterns a and /3 for a given randomness in 
Jxy ■ Assuming that the state a takes place with a probability W a for a given randomness, Parisi H i l27(i [ 
introduced a probability distribution of overlap values as 


P(q) = W a W^5(q - q ap ). 
a,(3 


(136) 


It is generally assumed that q aoi is independent of a in the limit TV —► oo. In that case, q aa = qEA and P(q) 
exhibits two Dirac-delta peaks at ±qEA, referr ing t o the global spin reversal symmetry of the Ising model. 
For the Sherrington-Kirkpatrick model, Parisi 27Cj 271 1 has shown the hierarchical (tree-like) structure in 
the overlapping q a p for a given randomness in J xy and discussed the general features of P(q) when averaging 
over the randomness in J xy . It was found that (( P(q)))j > 0 if —qEA < q < qEA- 

For more realistic spin glass models with short range interactions, basically different results have been 
reported by Newman and Stein 272 1. Middleton [273] . Billoire et al. 274 ] who studied the averaged overlap 
probability distribution. The differences are related to the number of energy valleys, the consequences of 
spatial structure and the absence/presence of frustration. At the same time these investigations raised many 
new questions. For a brief sur vey o f the differenc es an d recent achievements we suggest reading the recent 
reviews by Newman and Stein [275 ] and by Read [276 ], 

The effect of random-held h x on the the rmod ynamical behaviors of t he d -dimensional ferromagnetic Ising 
model was investigated by Imry and Ma [2771 ] , Villain 2781, Imbrie (279[ . using different approaches. In 
these investigations J xy = J, (h x ) = 0, and (h x ) = h? R are chosen where the mean variance h R measures the 
randomness. The complexity of questions related to the existence of ferromagnetic phase in the low noise 
limit is well reflected by the facts that the initial studies lead to controversial results. The systematic Monte 
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Carlo calculations (28fll. l28ll 


282l 283| have indicated the ferromagnetic phase in the presence of a weak 
randomness at low temperatures if d, > 2, otherwise there is no average magnetization in the presence of a 
random held. In other words, the paramagnetic to ferromagnetic critical phase transition can occur if hrt 
does not exceed a threshold value dependent on the lattice structures. Very recently Shrivastav et al. I284l | 
have investigated the morphological properties of the spin arrangement in the two- and three-dimensional 
systems in the low noise limit and reported power law behaviors in the correlation functions. In that case 
the paramagnetic phase possesses spin glass order qEA > 0, that increases when decreasing the temperature, 
although the system has only one paramagnetic state. 

Finishing this section we emphasize that the evolutionary potential games describing biological and 
social systems demand the systematic investigation of the effects of additional randomness that can occur 
in the connectivity structure, the local dynamical rules, and even in the number of strategies. On the other 
hand, the characteristic features of the spin glasses may be suppressed when slow variation is allowed in the 
parameters characterizing the randomness. 


8. ORDERING PROCESSES 

A typical two-dimensional ordering process in the Ising model from the random initial distribution 
towards the equilibrium state was illustrated previously in Fig. [15] This process can be observed for a fixed 
noise level below its critical value K c and it has some universal features. Namely, two equivalent types 
of ordered strategy arrangements form a domain structure that is topologically similar to the cases where 
islands are in lakes located within larger islands located in larger lakes etc. for the infinitely large systems. 
This pattern can be characterized by the average distance l(t) between two neighboring domain walls (along 
the horizontal or vertical cross-sections) that give a time dependent contribution to the equilibrium average 
potential. Namely, after a relaxation time, U eq — U(t) = a/l(t ) where a depends on the pair potential and 
noise level K. and l(t) oc \ft as long as l{t) < L. 

Figure [501 shows two typical Monte Carlo results obtained when the system with hawk-dove parameters 
is started from a random initial state on a square lattice and the time-dependence of the average potential 
[U(t)\ is recorded during the evolution for two noise levels. Despite the large system size (L = 4000) the 



Figure 30: Log-log plot of averaged U e q — U(t) vs. time for K = 0.1 (squares) and K = 0.34 (diamonds) at T = 1.5 and 
S = 0.2. The dashed line shows the slope (-1/2) characterizing the theoretical prediction. 

functions U{t) are decorated by fluctuations with an amplitude increasing with time. These undesired 
fluctuations are suppressed by averaging U(t) (at t = t* ~ 2*/ 2 ) over time intervals 0.8fi < t < 1.21, for the 
clear illustration of the asymptotic behavior mentioned. Notice, that it takes a longer time to achieve the 
asymptotic behavior when K approaches the critical point due to the critical slowing down as indicated by 
diamonds in Fig. [SH] 
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Different approaches have been developed to describe the above-mentioned time-dependent processes 
that remain valid also for the homogeneous three-dimensional systems. All the relevant features are well 
investigated by continuum description in terms of coarse-grained order parameter field [for det ailed discussion 
we suggest consulting the reviews by Hohenberg and Halperin [2851 ]. Langer [28fii ]. and Bray [287j ]. For this 
approach the systems are described by a scalar order parameter $ dependent on time and continuous space 
that can be considered as average magnetization (or strategy density) over a small region. 

Additionally, the motion of interfaces can be descri bed by geometrical appr oaches for both the two- and 
three -dimensional systems as it is detailed by Brakke |288l |, Brower et al. |289t |, and Goldstein and Petrich 
[29f)| (with further references therein). For these differential geometric descriptions the two-dimensional 
system is characterized by a set of closed curves representing the interfaces in the corresponding two-color 
map. The evolution of these curves is determined by a differential equation taking into consideration the 
average velocity of the interface depending on its local curvature, direction, and additional symmetries 
coming from the microscopic dynamical rules. For example, the reduction of the length of interfaces via the 
so-ca lled curvature driven interface evolution (for a very recent geometrical survey see the work by Garcke 
29l|) can explain different domain growing phenomena (l oc t a with a = 1/2, 1/3, or 1/4). Additionally, 


these approaches can describe a wide variety of interfacial instabilities occurring in solid state systems. On 
the other hand, the applicability of this method is limited to systems possessing one type of interfaces. 

In comparison to the spatial systems the final ordered strategy arrangement is formed significantly faster 
on small-world connectivity structures due to the absence of large distances. For the demonstration of this 
phenomenon we have performed MC simulations for an evolutionary hawk-dove game on a random bipartite 
regular graph with a degree of 4. The system started from a random initial state and the order parameter, 
the average payoff and potential were recorded after each MC steps performed at a noise level below its 
critical value. During the first steps the sharp increase of potential (see Fig. EU refers to the formation of 
ordered microdomains. Within this period the sublattice order parameter remains zero because of the large 



Figure 31: Average potential as a function of time during the sublattice ordering process on bipartite random regular graph 
for a hawk-dove game at T = 1.5, S = 0.5, K = 0.5, and N = 4 * 10 6 . MC data of four independent runs are illustrated with 
different types of lines. 

number of microdonrains of two types. The number of domains decreases gradually and after some time the 
size of the largest domains becomes comparable to the diameter of the given graph (oc In A''). Thereafter, 
the largest domain conquers the system. The parallel curves in Fig. [31] refer to similar scenarios and the 
randomness influences dominantly the time when ruling domain emerges. Evidently, on smaller systems 
the above behavior is disturbed by the stochastic events occurring both in the generation of the random 
networks and in the evolution of strategy distribution. 

8.1. Evolution in the limit K —> 0 

The pattern evolution in the Ising model on different lattices and graphs has been studied for several 
years in the zero noise limit |292t |. This particular case exhibits some curiosity. 


57 






For example, in the one-dimensional Ising model the long-range order can also be formed in the system 
if it is started from a random initial state while the evolution is controlled by the Glauber dynamics in the 
limit K — > 0 for h = 0. In that case the system becomes equivalent to the one-dimensional voter model 
29 -'ll 187t | and the interfaces (separating the opposite domains) move randomly. If two interfaces collide 


then both are annihilated. Finally the system evolves into one of the homogeneous states with a p roba bility 
dependent on the initial magnetization. The average domain size increases in time as l(t ) oc f 1 / 2 294 1. 

In general, the l(t) oc t 1 ! 2 scaling law is valid for the two- and three-dimensional lattices, too. However, 
in several trials the system evolves into a frozen state [295L ]296j and these events modify the long-time 
behavior. The bottom snapshots in Fig.[l2]show an example when rectangular boxes of the preferred phases 
are frozen into the opposite phase if the system is started form a random state being close to the ordered 
opposite structure. The application ofthe periodic boundary conditions may also result in frozen patterns 
in the finite systems. Spirin et al. cm have reported that for about one third of the trials in the two- 
dimensional lattice the system evolves into a poly-domain state where all the interfaces are horizontal (or 
vertical). For most of the cases only two strips are formed. It is reported furthermore that the number of 
frozen states is larger for d = 3. 

Recently Biswas and Sen 292| have studied the Ising system on a random network created from the 
one-dimensional lattice by adding new links into the connectivity structure. As a result the lattice sites 
have different degrees and for some constellations these irregularities were capable of blocking the domain 
growing processes, independently of the details of the generation of random netw orks. On the other hand, 
the freezing disappears in the limit N —> oo for densely connected networks 297 ]. 

The appearance of frozen patterns is expected for the n-strategy games (or Potts models) and also for 
systems where the number of neighbors is increased. 


8.2. Interfacial phenomena and rearrangement through nucleation 

Many relevant phenomena for the ordering or reordering processes can be well interpreted by considering 
microscopically the evolution of an interface separating two ordered phases for the two-dimensional ferro¬ 
magnetic Ising model. Figure [321 shows a horizontal interface between the up- and down-spin ordered phases 
that remains stable at K = 0 even in the presence of a weak magnetic field. Due to the rare stochastic 
events at low noise levels, one of the spins may reverse along the interface as denoted by the middle panel 
in Fig. 1321 The ferromagnetic nearest neighbor interactions enforce this spin to flip back, which represents 
the typical behavior. Sometimes, however, before the reconstruction of the smooth interface one of the 
neighboring spins may also reverse thus forming a two-spin cluster. 

The appearance of a two-spin cluster along this interface results in a new situation when the spins at 
sites x and y in Fig. [35] have similar environments, namely, there are two up-spins and two down-spins in 
their neighborhood and the preferred spin reversal is determined by the magnetic field. As a result this 
two-spin cluster can be considered as a nucleon from which the expansion (or growing) of the preferred state 
can start. The direction of the motion of steps and also its average velocity can be quantified by comparing 
the probabilities of spin reversals at sites x and y in the third panel of Fig. 1321 The quantitative analysis 
justifies that the average vertical velocity is proportional to the magnetic field (if \h\ <C 1). Throughout 



Figure 32: Consecutive steps in the evolution of a horizontal interface separating two ordered strategy arrangements on a 
square lattice. 


these consecutive elementary steps the system will evolve into the thermodynamically stable ferromagnetic 
state if 0 < K < K c . 
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The average vertical velocity of the above horizontal interface depends on two factors: the average 
horizontal velocity of the steps and the frequency of these steps along the interface. The second quantity 
depends on time and also on the nucleation rate characterizing the probability of the appearance of a 
sufficiently large nucleon. If the nucleation rate is low, then the interfaces are composed of large horizontal 
and vertical segments, otherwise the domain pattern is almost isotropic. 

Similar phenomena can be observed along the interfaces separating the equivalent anti-ferromagnetic 
domains. In that case the presence of the homogeneous external field h modifies the interface as it is 
illustrated in the snapshots of Figs. [T5] whereas the the average velocity of the interface is zero. At the 
same time, the application of a staggered magnetic field prefers one of the ordered spin arrangements to the 
opposite one and results in an average velocity proportional to h s (if h s -C 1 ). 

In the above process the appearance and expansion of the one-dimensional preferred nucleons play 
the crucial role in the evolution towards the final stationary state. Similar mechanism can be observed 
in d-dimensional lattices when the appearance of sufficiently large islands of the preferred ordered phase 
is necessary to initiate the transition from an ordered phase to the thermodynamically stable final spin 
arrangement. That happens when we wish to reverse the direction of magnetization by the application of 
an external magnetic field. In thermodynamical systems, the formation of the sufficiently large nucleon of 
the preferred phase is supported by a suitable series of stochastic elementary steps that may occur rarely. 
Sooner or later, however, some sufficiently large nucleons appear and catalyze the transition towards the 
stable phase via a domain growing process. The resident time in the meta-stable state(s) may be extremely 
long, especially at low noise levels. These phenomena are well investigated in a wide scale of physical systems 
and exploited in many products of high technology. 

In two-strategy evolutionary potential games the growth of the C domains on a square lattice can be 
characterized by the average velocity v of the moving step shown in Fig. 1321 If logit rule controls the 
evolution then 

e u*(C)/K e u v (D)/K 

V = e u x (C)/K + e u x (D)/K ~ e u y {C)/K + e u y (D)/K ‘ ( 137 ) 

Due to the similar strategy arrangements in the neighborhood the denominators are equal and 


g2(l +S)/K _ e 2 T)/K 
e 2(l +S)/K + e 2 T)/K 


( 138 ) 


in the social dilemma notation. As a result, within the region of stag hunt game the condition v = 0 defines 
a straight line (S = T — 1) on the T — S plane separating homogeneous cooperation (s x = C) and defection 
(sa, = D ) regions in the low noise limit (in agreement with the phase diagram plotted in Fig. 1161) . 


8.3. Interfacial phenomena in three- and n-state systems 

During evolutionary processes the interfacial phenomena can play crucial roles in the n-strategy (n > 2) 
systems, too. First we show what happens during the domain growth if the coordination type interaction 
between strategies 1 and 2 is extended by additional (neutral) strategies. Figure [ 55 ] illustrates the formation 
of two homogeneous domains (with strategy 1 and 2) if A = d(l,2) [see the definition (1411) at n = 5] if the 
two dimensional system is started from a random initial for logit rule at a noise level (K = 0.4). For this 
low noise level the strategies 3, 4, and 5 can occur very rarely inside the homogeneous domains. At the same 
time the mentioned strategies will be selected dominantly by the players located along the interfaces where 
the opposite effects of neighbors (with strategies 1 and 2 ) are balanced (see the sites x and y indicated in the 
right plot of Fig. [ 551 ) . In the present system the appearance of the additional strategies along the interfaces 
does not modify the general features of domain growing process. 

If the interactions are dominated by d(l,2) then the equivalence between the strategies 1 and 2 is 
evidently broken by the presence of self-dependent components and also by the additional coordination type 
interactions. More precisely, the preference of strategy 1 (or 2) depends on the value of 71 — 72 and also on 
the strengths of the components d(l, 3) and d(2,3). The comp etitio n between these components may result 
in different phase transitions as it is discussed by Vukov et al. j208| in a three-strategy evolutionary game. 
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Figure 33: (Color online) Two homogeneous domains of strategies 1 (white) and 2 (black) are growing in the two-dimensional 
system if the interaction is defined by A = d(l, 2) for n = 5. During the evolutionary process the strategies i = 3,4, and 5 
(denoted by red, blue and green boxes) occur dominantly at the steps of interfaces separating the white and black territories. 


For another example we mention the attractive n-state Potts model evolves towards one of the homoge¬ 
neous ordered states throughout a domain growing process at low noise levels, if the system is started from 
a random initial state. 

A similar phenomenon can also be observed for some spatial coordination games 31]. The left snapshot 
of Fig. [M] illustrates the spatial distribution of three strategies during the domain growing process in a 
system where A = d(l, 2) + d(2, 3) + d(l, 3) for n = 3]. The right hand snapshot of Fig. C2] shows similar 



Figure 34: (Color online) Typical strategy distributions on a square lattice for two types of the three-strategy coordination 
games during the domain growing process after 500 MCS if the system is started from a random initial state and the evolution 
is controlled by logit rule at low noise levels (K — 0.5A' C ). 

domain growing process in a system where the payoff matrix A' is obtained from A by exchanging its first 
and second rows. As discussed in Sec. cm and [311 ] there are two other equivalent systems that can be 
transformed into each other by exchanging two strategy labels in one of the sublattices. A similar symmetry 
is behind the equivalence of the ferromagnetic and anti-ferromagnetic Ising model on bipartite graphs. 

In two-dimensional spatial systems the geometrical features of these domain patterns differ significantly 
from those occurring in the two-state systems, where the interfaces form closed loops typically in the infinitely 
large systems (as mentioned in Sect. 0. For n = 3 one can distinguish three types of interfaces that can 
form closed loops or represent a planar network with three-edge vertices. 

In spite of the striking geometrical differences the average domain size increases with the square root of 
time (l(t) oc t 1 / 2 ) as reported by Grest et al. j298| and the geometrical features of these patterns become 
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similar o n th e scale of l(t) as it is observed for the Ising model and other field theoretical models surveyed 
by Bray 2871. 

In the present system, the existence of the third strategy does not influence the behavior of interfaces if 
the point defects are distributed sparsely. If the attractive Potts model includes a magnetic field favoring one 
of the homogeneous states then one can observe the expansion of the favored domains along their interfaces 
with an average velocity proportional to the magnetic field. 

A similar behavior is expected for three-strategy potential games if the payoff matrix g( 8 ) is weakly 
disturbed by additional self-dependent components (Idlil) . In those cases the difference between the corre¬ 
sponding 7 i values will determine the direction (and also the average velocity) of invasions between two 
’’homogeneous” domains. For example, this mechanism results in the prevalence of strategy 1 in the final 
stationary state if 71 > 72 , 73 . 

The right hand snapshot of Fig. [M] illustrates a system that is not yet studied systematically within the 
framework of Potts models. In the latter evolutionary game one can observe a homogeneous domain formed 
by strategy 3 and two equivalent chessboard like ordered strategy arrangements of the strategies 1 and 2. The 
most striking feature of this game is related to an inherent symmetry that we have discussed previously when 
justifying the equivalence of the ferromagnetic and anti-ferromagnetic Ising models on bipartite networks 
(see Sec. 17.41) . According to the generalization of the mentioned method, game g( 8 ) on a bipartite network 
becomes equivalent to those where the pair interactions are defined by g( 6 ) if the players in sublattice Y 
exchange the labels of their first and second strategies (1 -o- 2). The latter transformation is equivalent to the 
exchange of the first and second row of the payoff matrix. Due to the mentioned relation the three domains 
are equivalent on the square lattice and after a domain growing process one of these ordered structures will 
prevail in the finite systems. When increasing the noise level I\ this system undergoes an order-disorder 
critical phase transition belonging to the universality class of the three-state Potts model. 

One can generate two additional potential games with payoff matrices obtained from g( 8 ) by exchanging 
the labels 1 -O- 3 or 2 -o- 3 for the players staying in sublattice Y. The resulting games are similar to g( 6 ). 
These relatives of the Potts model can be constructed as suitable linear combinations of three elementary 
games (g( 6 ), g(7), and g( 8 ) j3lj. It turned out that the corresponding three-dimensional subset of games 
can be considered as a generalization of the Potts model. Here we have to emphasize that the additional 
self-dependent components are not capable of preferring one of the sublattice ordered two-strategy structures 
to its anti-pattern. 

The exploration of the three-strategy symmetric potential games is not complete. Even more complex 
behavior is expected for n > 3 when the number of interfaces as well as the types of vertices increases with 
n. The preliminary results h ave indicated the dominance of three-edge vertices that follows a complicated 
transition/annihilation rule l299l |. The consideration of nonsymmetric games will allow us to study the effect 
of those types of self-dependent components that distinguish the chessboard and anti-chessboard arrange¬ 
ments of two strategies, as it is done when applying a staggered magnetic field to the anti-ferromagnetic 
Ising model. 

In the above-discussed models the average motion of the invasion fronts is driven by the increase of 
individual payoffs that is quantified by the increase of U in the thermodynamic potential $ (15121) if K —> 0. 
In the maximization of d>, however, the high entropy of the disordered phase, especially for a large number 
of strategies, can become the leading term at sufficiently high values of K. 

Figure l35l illustrates a domain growing process on a square lattice from one of the ordered phases (here 
s x = 1) to the disordered strategy arrangement if the interaction is described by d( 1, 2) for n = 50 at a noise 
level K = 0.52 > iv c (50) = 0.512(2). The process begins with a nucleation procedure that is followed by the 
expansion of the disordered territories where all the strategies are present with approximately equivalent 
probabilities. 

Here the evolution is also controlled by the logit rule, therefore we can estimate the average velocity v 
of a step along the interfaces in the same way as it is described in Sec. 18.21 Neglecting the appearance of 
strategy 2 the value of v can be approximated as 


o2/K 


n — 2 


v ~ 


z 1 ! K + n — 2 e 2 / K + n — 2 


(139) 
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Figure 35: (Color online) Three consecutive snapshots at times t = 200, 1000, and 3000 MCSs (from left to right) show a 
domain growing process when the (white) homogeneous spatial strategy distribution (s^ = 1) transforms into the disordered 
phase composed of n = 50 strategies (distinguished by different gray scales). 


Accordingly, the disordered phase expands at the expense of ordered phase if v < 0, that occurs if K > 
2 / In (n — 2) for sufficiently large values of n when the /v-dependence of strategy frequencies exhibits a first- 
order phase transition. The criterion v = 0 gives us an estimation for the critical point, K c (n) = 2/ In (n — 2), 
that agrees very well with the Monte Carlo result given above for n = 50. 

The latter result implies that the high-entropy phases can occur for any noise level if n is sufficiently 
large. Here it is worth mentioning that similar arguments justify the stability of high-entropy alloys at room 
temperatures. Recently the high-entropy alloys are studied pr ogres sively and considered to be a promising 
family of materials with a wide scale of applications (300l . l30ll l302j . 

Anyway, the average velocity of an interface can be determined numerically if the Monte Carlo simulations 
are started from artificial initial state where the competing phases are represented by two domains with equal 
sizes . Us ing this method one can evaluate the phase boundaries in the phase diag rams more accurately 
l.'SO.l 1304 ], particularly if the first order transition is accompanied with a hysteresis ]305l] or sensitivity to 
the initial state [3061] . 


8-4- Slow relaxation in random systems 

Up to now we have mainly discussed the stationary states of the systems (in the limit N —> oo). In 
general, homogeneous systems evolve towards the stationary states exponentially if the state is weakly 
perturbed. We have mentioned two exceptions when the homogeneous system reaches the final state more 
slowly. In the first case the system evolves form a random initial state into one of the ordered arrangements 
through a domain growing process and the average domain size (or correlation length) increases with \/t as 
detailed above. The second case occurs at the critical point where system behavior is dominantly controlled 
by the fluctuations and results in a power law decay in most of the quantities. Now we briefly discuss the 
slow relaxation processes observed in the Griffiths phase of the random Ising systems. 

The slow (nonanalytic) relaxation of the magnet izati on in the paramagnetic phase of a random ferromag¬ 
netic Ising model at h = 0 was reported by Griffiths [3071 who studied a diluted Ising model on a lattice where 
a portion of the lattice sites are not occupied by Ising spins id08 . Subsequent analyses have indicated the 
presence of Griffiths phase in many other random Ising models between the paramagnetic and ferromagnetic 
(or spin-glass) phases. According to the investigation of different models the one-site two-time correlation 
function [defined in Sec. [8] is found to have a ”stretched-exponential” form, c/(0,t ) ~ exp| — (f/r) K ] with 
0 < k < 1, depending on the spatial dimension and other details of the system 30$ . 3ld 311 . 312 ], 

The related spatial patterns assume the existence of sparse and large ordered domains in the random 
environment. The tliermalization (spontaneous reversal) of these sparse domains is very slow because it 
requires a long sequence of coherent flipping over a large volume. According to this picture the time- 
dependence of magnetization (or any order parameter) can be approximated as 


/•OO 

m(t) = / w{t) exp[— t/r\ dr 

Jo 


(140) 
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with a suitable choice of the weight function w(t). Noest [313L l314l| has shown that if the relaxation time 
increases exponentially with the size n of a compact cluster (r ~ exp [an]) and the probability of such clusters 
decreases exponentially with n, then the leading term of the asymptotic behavior of rri{i) can be estimated 
as 


m(t) ~ t 9 


(141) 


where 9 (6 > 0) depends on a and on other parameters within the Griffiths phase. It is emphasized that 
similar behavior is reported for several other system s, e.q ., stochastic cellular autom ata |314l| and contact 
processes with quenched disorder in the environment |315| . Today the contact process [31f| is considered as 
the paradigm of systems where the extinction of a species/strategy exhibits a critical transition that belongs 
to the directed percolation universality class .317 . For a survey of the main features of this critical transition 


318 . 


In the latter system the quenched randomness modifies 
319|. The presence of Griffiths phase and its consequences 


we suggest consulting the review by Hinrichsen 
also the system behavior at t he cr itical point 

are described by Munoz et al. [320] for the contact process on complex networks. 

The occurrence of the Griffiths phase in the contact process has implied algebraic extinction of a strategy 
in many evolutiona ry g ames where the evolutionary rule is based on imitation in spatial systems with 


quenched disorder [32l|. In most of the evolutionary games with quenched randomness, however, the 


appearance of Griffiths phase is not investigated although it can cause incorrect numerical data in the 
vicinity of the critical point (s). 

The above theoretical picture supposes that the relaxations of the domains are independent of each other 
(i.e., the rare events are not organized hierarchically). This feature simplifies the numerical analysis of these 
systems as the Monte Carlo simulations can be performed simultaneously on many ’’sm all” systems. In the 
opposite cases, when the models involve hierarchically constrained dynamical processes [.309t |. more complex 
finite-size analyses are required. 

At the end of this section we underline the relevance of the Griffiths phase in the evolutionary games 
modeling biological or social systems where the quenched randomness is assumed naturally [[260]. For the 
numerical analysis of these systems in the Griffiths phase we have no chance to achieve the final stationary 
state due the the slow algebraic relaxation processes. It means, on the one hand, that in the Griffiths phase 
the final stationary quantities should be determined by extrapolation of the asymptotic behavior. On the 
other hand, for the interpretation of the experimental and numerical data we should consider the fact that 
the system has not achieved its equilibrium state. 

The Griffiths phase represents technical difficulties in numerical simulations that can be avoided by 
introducing slow variation in the randomness that is also a natural ingredient of biological and social systems. 
The relevant differences between the quenched and temporal randomness are detailed by [318j for the contact 
process. 


9. DEVIATIONS FROM THERMODYNAMICAL EQUILIBRIUM 

In the absence of potential the existence of Boltzmann distribution becomes meaningless and we cannot 
apply the results of equilibrium statistical physics and thermodynamics. Additionally, the validity of ther¬ 
modynamics is dropped when applying an evolutionary rule that breaks the detailed balance and drives the 
system far away from the Boltzmann distribution, as happens, for example, in the imitation-based dynamics 
used frequently in many previous investigations. For some types of coevolutionary games the absence of 
the fixed connectivity structure, the possible changes in personal features and spatial location raise many 
additional questions whose analyses go beyond the scope of the present work. 

In the following sections we discuss briefly some effects of the matching pennies and rock-paper-scissors 
games that can be studied in multi-agent evolutionary games even for the application of logit rules. The 
discussion of these games is unavoidable because of their important role in the subgames of any nxn matrix 
games, which can affect the system’s behavior significantly. 

9.1. Effects of matching pennies 

The matching pennies game, defined by f'(8) in its bimatrix form (15H1) . represents the simplest cyclic 
interaction. For a two-player evolutionary game with logit rule the unsatisfied player reverses her strategy 
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with a high probability dependent on K , while the opposite transitions become rare. Thus this interaction 
breaks the detailed balance and can be considered as a microscopic force inducing cyclic variations in 
the systems. The effect of this interaction can be quantified by the probability current(s), measuring the 
difference in the frequency of forward and backward transitions along the directed edges of the flow graph 
(see Fig. [5]). This current is uniform along the four edges in the stationary state of this stochastic process. 

The induced circular probability current creates observable variations in the probability of the strategy 
profiles. For the illustration of this phenomenon we discuss a system where a weak matching pennies 
component is added to the payoff of a hawk-dove game. The variations in the probabilities are illustrated 
in Fig. [36] In this figure the height of the dark columns denotes the probability of each strategy profile in 



Figure 36: Probabilities of the four strategy profiles are proportional to the height of the dark (e = 0) and white (e = 0.1) 
columns for a two-person hawk-dove game if it is extended by a matching pennies component with a strength of e. The arrowed 
gray circle denotes the probability current loop induced by the matching pennies game. 


the stationary state when G = G^ sd ^ as defined by (1551) . For the given parameters ( T = 1.4, S = 0.3, and 
K = 0.3) p(l, 1) < p{ 2, 2) whereas p(l, 2) = p( 2,1). If the latter game is modified by adding a weak matching 
pennies component to the payoffs (quantitatively G = G( sd ) +£G^ mp ^ with £ = 0.1) then the appearance of 
probability current is accompanied with a striking variation in the stationary state. The largest increase of 
p{ 2,1) can be interpreted as the consequence of a congestion phenomenon. Accordingly, for the maintenance 
of the circular probability current, the lowest value of p(l,l) plays the role of a narrower bottleneck that 
creates an increase in p{ 2,1) to ensure uniform probability current through the four-edge loop. At the same 
time the value of p( 1, 2) is decreased, that is, the presence of the matching pennies component destroys the 
equivalence of the two Nash equilibria. 

In the multi-agent version of this evolutionary game on the square lattice, the above microscopic effect 
occurs for each interacting neighbor and affects the macroscopic behavior. This extension of the models can 
be performed for those lattices that can be divided into two sublattices (X and Y, where the two types of 
players are located separately). As a result of this breaking of the original symmetry, one of the sublattice 
ordered strategy arrangements can be preferred in the low noise limit. Evidently, the preference is reversed 
with the sign of e. Furthermore, the preference is also reversed if we consider the upper half part of the 
hawk-dove game, where S > (T — 1) and p{ 2,2) < p(l, 1), because here the strategy pair (2,2) plays the 
role of the narrower bottleneck in the circulation. In fact, the effect of the matching pennies component is 
similar to the application of a suitable staggered magnetic field in the anti-ferromagnetic Ising model. 

Along the line S = (T— 1) in the parameter space, the above effect does not work because here p{ 1,1) = 
p( 2,2). In these systems the presence of th e matching pennies component does not destroy the universal 
features of Ising type critical transition [3221 ], whereas the value of K c is reduced (proportionally to e) if its 
strength does not exceed a threshold value. In the latter case the pair interaction belongs to the classes of 
ordinal potential games, because weak contribution of the matching pennies is not enough to change the 
edge directions in the flow graph. 

In non-equilibrium systems the breaking of detailed balance can be well quantified by evaluating the 
entropy production / (for a survey see [411] b This quantity is constructed from the frequencies of the 
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forward and backward transitions between states s and s' if these transitions are allowed by an elementary 
step in both directions. Now the random consecutive elementary steps consist of single site strategy changes 
from (s x ,s_ x ) to (s x ,s_ x ) appearing with a frequency W{s x —► s x ) in the stationary state. The entropy 
production summarizes the contributions of the forward and backward transitions for each possible transition 
pair in the following way: 

1 = 2 ^2 l w ( s * ~ w ( s z s ^] ln irls^Tsty' ( 142 ) 


Notice that this quantity is always positive (/ > 0) except the case 1 = 0 when the conditions of detailed 
balance are satisfied, i.e., if W{s x — > s' x ) = W{s' x — > s x ) Vs x , s' x , s_ x . Despite the large number of transitions 
in lattice systems the entropy production can be well estimated by recognizing that the transition frequency 
IF(s x —> s' x ) depends dominantly on the close neighborhood of player x, who modifies her strategy, if the 
dynamics is controlled by short range interactions. Besides it we can exploit the translation invariance of 
the lattice system. As a result the specific entropy production (/ /N) can be estimated by considering the 
transition frequencies at any site x for all possible strategy configurations in its close neighborhood. 

For example, during the Monte Carlo simulations on a square lattice one can determine the transition 
frequencies s x —> s x for all the 2 4 = 16 possible strategy configurations of the four nearest neighbor sites and 
also for the cases when 2 8 = 256 configurations are distinguished on the first- and second neighbor sites. In 
this way we can deduce two approximate results for the specific entropy production ( I/N ) and comparison 
of them indicates the relevance of the second neighbors although they do not influence directly the transition 
in the present models. Evidently, the larger the neighborhood, the more accurate is the present approach 
(for a more detailed description of this approach see @). 

Figure [37] shows the Monte Carlo results for the specific entropy production when varying the noise level 
for different strengths of the matching pennies component. Here the closed symbols represent data obtained 
when only the first neighbors are taken into consideration in the identification of s_ x . Notice that the latter 
approximate data are close to those we obtained when a larger neighborhood (the first and second neighbors 
of player x) is used for the characterization of s_ x . The small differences between the two sets of data 
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Figure 37: Specific entropy production as a function of K for evolutionary game with G = Gl sd l + pair interactions 

on the square lattice if T = 1.4, S = 0.3, and e = 0.1 (open diamonds and bullets), 0.175 (open boxes and closed diamonds), 
0.25 (open circles and closed boxes). The open symbols denote data obtained for the larger neighborhood. 

justify the reliability of this approach. Notice furthermore another general feature: I/N vanishes in the 
limit K —> oo when only the randomness controls the players’ decisions. 

Figure [37] illustrates that the force of ordering (strength of hawk-dove component) blocks the strategy 
reversals and also the breaking of detailed balance in the low noise limit if e is less than a threshold value, 
more quantitatively, |e| < \eth\ where \eth\ = l/2min(|T— 1|, liS)). In the opposite limit, the specific entropy 
production diverges if K —> 0 and we can observe a remarkably different behavior. 
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The divergency of I/N in the low noise limit is characteristic to systems where cyclic dominance controls 
the system behavior and prevents the formation of ordered strategy arrangements. A similar divergency 
occurs in systems where a finite portion of the transition pairs becomes unidirectional. If the evolution is 
defined by only the matching pennies component in a spatial evolutionary game then the visual observation 
indicates random strategy distribution on a square lattice. The qua ntita tive analysis, however, has clearly 
indicated weak correlations between the second and third neighbors [322j | that evidently vanish in the high 
noise limit. 


9.2. Effects of rock-paper-scissors game 

If the multi-agent evolutionary games are composed of equivalent players with three strategies, then the 
rock-paper-scissors component is responsible for the deviation from the potential games. The rock-paper- 
scissors game itself creates a weakly correlated random distribution of the three strategies on a square lattice 
if a logit rule controls the evolution. The three strategies are present with the same probability (1/3) and the 
numerical investigations indicate a weak spatial correlation that is similar to those found for the matching 
pennies game. The latter analogy implies that the alternation of the three strategies on each site reflects 
relevant breaking of the detailed balance, particularly in the limit K —> 0 when I/N diverges. 

Contrary to the two-strategy games with a matching pennies component, the presence of a weak rock- 
paper-scissors component can cause more relevant changes in the macroscopic behavior of the multi-agent 
spatial three-strategy evolutionary games. To demonstrate this, we consider the three-state Potts model 
on a square lattice if the uniform attractive pair interactions are modified by introducing a weak cyclic 
dominance, i.e., G = G( Potts ) + eG^ where G( Potts ) = d(l,2) + d(2,3) + d(l,3). As mentioned before, 
the n = 3 Potts model evolves into one of the three ordered states at low noises. In the presence of a 
weak cyclic dominance, however, the domain growth is blocked by the formation of rotating spiral arms 
as illustrated in Fig. [35] The spiral form of the rotating edges is a direct consequence of the fact that 
the component eG^ rsp ) induces an average invasion velocity independent of the direction of the interface 
and of the distance measured along the interface from the center of the given three-edge vortex (vortex 
means rotating vertex). Due to the cyclic symmetry, all the three edges of a vortex (anti-vortex) rotate 
clockwise (anti-clockwise). Sometimes the moving interfaces meet and may annihilate each other or create a 
new vortex-anti-vortex pair. The latter processes are accompanied with a rearrangement of the connections 
(represented by the interfaces) between the vortices and anti-vortices. Some geometrical features, namely, 
the average curvature of interfaces, the average distance of vortices, the average le ngth of interfaces bet ween 
a vortex and anti-vortex pair, were already investigated by Szolnoki and Szabo [323J, Szolnoki et al. 324 


who found that the critical transition is suppressed in the presence of cyclic dominance and the correlation 
length diverges as oc 1/e if e —> 0 at sufficiently low noise levels. 



Figure 38: (Color online) Rotating spiral arms in the snapshots characteristic to the self-organizing pattern on the square 
lattice where the pattern evolution is controlled by attractive (ferromagnetic) three-state Potts type interaction and a rock- 
paper-scissors component with a strength of e = 0.1 for logit dynamics at a low noise level. 
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If the previous model is modified by additional potential game components that favor one of the homo¬ 
geneous states, then a sufficiently weak cyclic component cannot prevent the formation of an ordered state 
at low noises. In these cases the interaction itself might be considered as an ordinal potential game. A 
further increase in the cyclic component, however, can maintain a self-organizing spatio-temporal pattern 
(see Fig. 1551) in which the portions of the three strategies are different. In the practice of evolutionary game 
theory this mechanism is exploited when the undesired effe ct of social dilemma is reduced by introducing a 
third strategy representing voluntarism [325[, ” tit for tat” 326 ], or punishment [327 1. 

It is emphasized t hat simila r rotating spirals are observed in ma ny o ther systems including 
Zhabotinsky rea ction [328l l329l| . excitable media ( e.g. cardia c mu scle [33(1 . neural systems l.i.'l 11 ) 


Belousov- 
_ _ epidemi¬ 

ological models 13321 . and biological/ecological models 333. 133 ll 335|. The robustness of similar spatio- 
temporal patterns is also demonstrated by numerous three-strategy spatial evolutionary games (for a survey 
see [22j). It is already well-known that the rock-paper-scissors type cyclic dominance helps the maintenance 
of all the participating strate gies/speci es 162 even for inhomogeneous invasion rates generating a non-trivial 
reaction in the populations [33fil 13371 ]. The cyclic dominance can mediate positive or negative feedback 
thro ugho ut the cyclic process that depends on the parity of the number of strategies within the cyclic game 
338l 339l . 340[ . The consequences of this parity effect can affect the behavior in many systems. 


The presence of cyclic interactions is re sponsible for the survival of a large number of strategies and 


the bio-diversity in biological systems [3411 . 22]. The investigations of predator-prey models with a large 


number of species and with a complex structure of cyclic dominance indicate an extremely wide scale 
of behaviors. In some systems the cyclic components can maintain different subsets of strategies, called 
strategy associations, that are stabilized against the external invaders by suitable spatio-temporal patterns. 
These strategy ass ociations can survive simultaneously in a large spatial systems by forming large domains 
342l 343l . 344 . 1345 1. In these complex systems the cyclic dominance among t he st r ateg y associations can be 
quantified by determining the average velocity of interfaces separating them [3461 . 13471 ]. 

The analysis of the competition between strategy associations is generally based on models with random 
sequential imitation type evolutionary rule that may even be applie d sim ultaneously for several systems |348| 
when the models become similar to stochastic cell ular automata [3491 . l350j |. The spatial rock-paper-scissor 
game with a synchronized stochastic logit update 135 ll ] has demonstrated the app earance o f chim era state s 
which have been intensively studied in the literature of coupled spatial oscillators [352 , 353 . 354 3551 356 1. 
For the repeated two-player rock-paper-scissors game the synchronized logit rule at low noises results in 
cyclic choices [e.g. (1,1) —> (2, 2) —> (3, 3) —> (1,1)] until the first mistake. Afterwards the (1, 2) —> (3, 2) —> 
(3,1) —> (2,1) —> (2,3) —> (1, 3) —> (1, 2) cycle is repeated in the absence of mistakes. Such cycles can also 
occur in the spatial system as it is illustrated in Fig. 1391 The numerical simulations indicate the formation 




Figure 39: (Color online) Snapshot (left) of a chimera state on a square lattice where the evolution of a rock-paper-scissors 
game is controlled by a low noise logit rule applied synchronously. The right hand plot illustrates how the spatial patterns 
change cyclically. 


of large domains in which the cyclic choices are stabilized by the neighborhoods suppressing the effect 
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of individual mistakes. In the snapshots one can distinguish nine types of domains representing different 
phases of these cycles. Additionally the visualization of the pattern evolution indicates clearly the presence 
of rotating spiral arms due to the cyclic dominance between the oscillating associations. 

Finally we mention that the interfaces separating two strategy associations can serve as a location for 
the emergence of a new strategy association with a proper spatio-temporal structure. Such phenomena were 
reported in spatial evolutionary games with cyclic dominance between the strategies/spe cies (for n < 5) 
where the evolutionary rule is based on imitation and site exchange of neighboring players 1 51. 'll |346| . The 
interfaces in Fig. [33] can also exemplify the appearance of a new phase that may play a crucial role in the 
formation of multicellular living materials. 

10. CONCLUSIONS AND OUTLOOK 

We have reviewed our recent understanding of potential games representing an intimate relationship be¬ 
tween physical systems and models applicable to study relevant phenomena in biological and social/economical 
systems. The analogy becomes particularly important for those social and biological multi-agent systems 
where the pair interactions can be well described by symmetric n x n potential games with logit rules when 
the systems are driven into the Boltzmann distribution and the general laws of thermodynamics are valid. 
The application of concepts and methods developed in statistical physics proved to be trivially beneficial 
for the partnership games where the equivalent players share their income equally. It turned out, however, 
that the analogies can be directly extended to the potential games representing a wider scale of games. 

The evaluation of the potential is demonstrated if it exists. The largest component of the potential matrix 
identifies the preferred Nash equilibrium (playing the role of ground state in physical systems) even for multi¬ 
agent systems composed of uniform pair interactions. This feature implies a simple method for determining 
the phase diagram at low noises. If the largest component of the potential matrix is located on the main 
diagonal of the potential matrix then all the players prefer to choose the corresponding strategy independent 
of the connectivity structure. In general, the latter systems show a thermodynamical behavior represented 
by the Ising model in the presence of an external magnetic field, even for n > 2. As Ising type models have 
already been investigated under a wide scale of different conditions (including symmetries, randomness, 
networks, ordering processes) therefore many results of statistical physics can be directly adapted to explain 
the phenomena in evolutionary games, as well. If the largest pair of the potential matrix components occurs 
outside the main diagonal then the systems have two equivalent preferred Nash equilibria and become 
similar to the anti-ferromagnetic Ising models. On bipartite graphs these systems exhibit an Ising type 
order-disorder phase transition when the noise level is increased, otherwise the sublattice ordering can be 
suppressed by frustration and/or randomness that may result in extremely slow relaxations towards the final 
stationary state as discussed in the literature of spin glasses and Griffiths phases. 

The classification of games into four classes of interactions has helped us conclude general features 
characterizing the corresponding subset of games. Accordingly, the players are not interested in favoring 
one of their strategies for games with cross-dependent payoffs, thus they choose strategies at random. For the 
self-dependent payoffs the players’ behaviors can be considered separately from each other. Consequently, 
all these multi-agent systems with equivalent players and interactions are well described by considering only 
one player. It is found, furthermore, that the real pair interactions of potential games can be built up as a 
linear combination of coordination type games between the possible strategy pairs. 

It is shown that the presence of the elementary games with cyclic dominance prevents the existence 
of potential. The consequences of the latter deviations are discussed for several examples challenging the 
application of methods developed in the field of non-equilibrium statistical physics. In the light of these 
results one can conclude that these terms of interactions result in self-organizations characterizing living 
systems. 

We have shown that some general questions of traditional game theory become particularly transparent 
when using the tools of graph theory. Here we used graphs for three purposes. The dynamical graphs 
visualize the strategy profiles (microscopic states) and the possible transitions between them if only unilateral 
strategy changes are allowed in the system. Due to the simple structure of the dynamical graph (for the 
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symmetric games) we could determine the number of independent and relevant loops, along which the sum 
of payoff variations should be zero for the potential games. The flow graph illustrates the preferred strategy 
changes and simplifies the identification of the pure Nash equilibria (that always exist in potential games) as 
nodes without outgoing edges. The dominance graph denotes graphically the payoff differences quantified by 
the antisymmetric part of the payoff matrix. Using these concepts one can distinguish cyclic and hierarchical 
dominance. In potential games the hierarchical dominance can be related to emergence of social dilemmas 
occurring even for n > 2. It is hoped that further graph theoretical investigations can throw light on 
additional relationships. 

When writing this survey we faced challenging questions and interesting phenomena week by week. Some 
of these problems have already been clarified during the preparation of this work while others remained in 
the state of ” challenging questions”. Examples for the former problems are the decomposition of matrix 
games, the identification of different classes of interactions, the relevance of cyclic dominance, the inherent 
symmetries involved in the matrices, and a series of interesting phenomena. The list of the latter examples 
is much longer and contains the identification of ordinal potential games, the elucidation of inherent symme¬ 
tries in the classes of interactions, the systematic investigation of social dilemmas for potential games, the 
relevance of high-entropy associations, the spontaneous formation of strategy associations in the presence 
of cyclic games, the extension of relaxation process for quenched random interactions, the co-evolutionary 
processes including the evolution of connectivity networks, payoffs, dynamical rules and emergence of new 
strategies, etc. The study of these intriguing questions offers further promising challenges. 
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